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This  document  represents  a final  report  to  the  Office  of  Naval 
Research  on  Task  No.  MR  064-554,  the  Display  of  Complex  Three  Dimen- 

" 1 ■!  II  ,||  | . 

sional  Finite  Element  Models  - Phase  III.  The  main  body  of  the  report 
is  a discussion  of  the  computer  programs  which  have  been  developed  dur- 
ing Phase  III.  This  material  was  prepared  by  Thomas  W.  Sederberg  and 
serves  as  his  M.  S.  Thesis  at  Brigham  Young  University.  The  remain- 
ing paragraphs  in  this  forward  discuss  activities  related  to  the  con- 
tract that  are  not  discussed  in  the  main  body  of  ;!ie  report. 

Distribution  of  MOVIE. BYU 

This  general  purpose  computer  graphics  software  package  (largely 
the  result  of  efforts  under  Phases  I and  II  of  this  contract),  has  now 
been  distributed  to  approximately  160  organizations  in  the  United  States, 
Canada,  England,  France,  Germany,  Norway,  Israel,  and  Australia.  A 
complete  mailing  list  of  those  organizations  is  included  as  Appendix  C. 

Technical  Presentations 

Since  the  award  of  Phase  III  of  this  contract,  technical  pres- 
entations featuring  results  obtained  under  this  funding  have  been  made 
by  Dr.  Christiansen  to  the  following  groups.  Of  cQucse,  travel  funds 
came  from  many  sources  and  some  of  the  presentations  were  on  an  informal 
basis. 

1.  Envirotech  Corporation  (Winston-Salem,  North  Carolina) 

2.  Waterways  Experiment  Station,  Corps  of  Engineers  (Vicksburg,  Miss.) 

i 

3.  University  and  College  Designers  Association  - National  Conference 

(Provo,  Utah)  ! 
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4.  Engineering  Society  of  Detroit's  Third  Annual  Computer  Graphics 


Conference  (Detroit,  Michigan) 

5.  Industrial  Design  Department,  Center  for  Creative  Studies  (Detroit, 
Michigan) 

6.  Computer  Graphics  Workshop,  University  of  Arizona  (Tucson,  Arizona) 

7.  Graduate  Seminar  in  Computer  Science,  Brigham  Young  University 

(Provo,  Utah)  ' 

8.  3rd  Southwest  Graduate  Research  Conference  in  Applied  Mechanics, 
University  of  Texas  (Austin,  Texas) 

9.  Department  of  Creative  Arts,  Purdue  University  (Lafayette,  Indiana) 

10.  Raytheon  Missile  Systems  (Bedford,  Mass.) 

11.  Graphics  Utah  Style  - 77  (Snowbird,  Utah) 

12.  Symposium  on  Computer  Methods  in  Engineering,  University  of  Southern 
California  (Los  Angeles,  California) 

13.  Genisco  Corporation  (Irvine,  California) 

14.  Tektronix  Corporation  (Wil sonvil le,  Oregon) 

15.  Shell  Development  Co.  (Houston,  Texas) 

16.  Design  and  Drafting  Seminar,  Brigham  Young  University  (Provo,  Utah) 

17.  American  Society  of  Mechanical  Engineers,  Annual  Meeting  (Atlanta, 
Georgia) 

18.  Art  Directors  Club  (Salt  Lake  City,  Utah) 

Since  the  final  report  for  Phase  II,  the  following  technical  papers 
have  been  published. 

Christiansen,  H.  N.,  Brown,  B.  E.,  and  McCleary,  L.  E.,  "A  General  Purpose 
Computer  Graphics  Display  System  for  Finite  Element  Models,"  46th  Shock 
and  Vibration  Bulletin,  Part  5,  August  1976,  pp.  61-66. 
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Christiansen,  H.  N.,  "Computer  Graphics  - Treatment  for  the  Terminal 
Illness,"  Preprint  2765  - Development  of  Computational  Methods  in  Str- 
uctural Analysis  and  Design:  Past,  Present,  and  Future,  ASCE,  Philadelphia, 
PA,  Sept.  1976,  pp.  VI  1-11. 
t Christiansen,  H.  N.,  "Computer  Simulation  of  Distorted  Structural 

Frameworks,"  Journal  of  Computers  and  Structures,  Vol . 6,  Dec.  1976, 

pp.  497-501. 
t 

f Christiansen,  H.  N.,  and  Stephenson,  M.  B.,  "MOVIE. BYU  - A General 

Purpose  Computer  Graphics  Display  System,"  Proceedings  of  the  Symposium 
on  Applications  of  Computer  Methods  in  Engineering,  USC,  Aug.  1977. 

Recently,  three  more  technical  papers  (two  covering  aspects  of 
Phases  I and  II  and  other  concerned  with  Phase  III)  have  been  written 
and  submitted  for  upcoming  meetings. 
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Chapter  1 

INTRODUCTION 

There  has  been  a disparity  between  the  conventional  method 
of  describing  topographic  surfaces  (i.e.  contour  line  definition) 
and  a format  of  surface  description  often  used  in  continuous-line 
computer  graphics  (i.e.  panel  definition).  The  two  differ  enough 
that  conversion  from  contours  to  panels  is  not  a trivial  problem. 

A computer  program  that  performs  suci  a conversion  would  greatly 
facilitate  continous  tone  display  of  topographical  surfaces,  or 
any  other  surface  which  is  defined  by  contour  lines. 

This  problem  has  been  addressed  by  Keppel'  and  alluded  to  by 

2 

Fuchs  . Keppel's  is  an  highly  systematic  approach  in  which  he  uses 
graph  theory  to  find  the  panel  arrangement  which  maximizes  the  volume 
enclosed  by  concave  surfaces.  Fuchs  mentions  an  approach  to  the 
problem  as  part  of  an  algorithm  to  reconstruct  a surface  from  data 
retrieved  from  a laser  scan  sensor. 

This  thesis  elaborates  on  a general  conversion  system. 
Following  a brief  overview  of  computer  graphics,  a simple  algoritm 


E.  Keppel , "Approximating  Complex  Surfaces  by  Triangulation 
of  Contour  Lines,"  Journal  of  Research  and  Development,  IBM  Vol . 19, 
No.  1 (January  1 975) > 2-11. 

2 

Henry  Fuchs,  "The  Automatic  Sensing  of  3-Dimensional  Surface 
Points  from  Visual  Scenes"  (unpublished  PhD  dissertation,  University 
of  Utah,  1975.) 
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, is  described  which  extracts  a panel  definition  from  a pair  of 

* 

adjacent  contour  loops  subject  to  the  restriction  that  the  two  loops 
are  similarly  sized  and  shaped,  and  are  mutually  centered.  Next,  a 
mapping  procedure  is  described  which  greatly  relaxes  the  above 
restrictions.  It  is  also  shown  that  the  conversion  from  contours  to 
panels  is  inherently  ambiguous  (to  various  degrees)  and  that  occa- 
sionally the  amibiguity  is  great  enough  to  require  user  interaction 

% 

to  guide  the  conversion  algorithm.  An  important  complication  add- 
ressed in  this  thesis  is  the  problem  of  handling  cases  where  one 
contour  loop  branches  into  two  or  more  (or  vice  versa). 

Attention  turns  next  to  a contour  line  definition  of  the  human 
brain,  and  special  problems  encountered  in  preparing  those  data  for 
continuous  tone  display.  The  final  chapters  explain  the  fortran 
implementation,  present  an  example  problem,  and  show  sample  pictures 
of  the  brain  parts. 
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Chapter  2 

AN  OVERVIEW  OF  COMPUTER  GRAPHICS 

The  past  decade  has  seen  fantastic  advances  in  the  field  of 
computer  graphics.  Today,  it  is  a sheltered  person  who  is  not 
familiar  with  some  form  of  computer  graphics,  be  it  Snoopy  calendars 
or  computer  ping-pong  on  one  end  of  the  spectrum,  or  sophisticated 
airline  pilot  training  simulators  on  the  other  end.  Display  mediums 
used  in  graphics  are  very  diverse,  and  include  raster  scan  cathode 
ray  tubes,  cathode  ray  storage  tubes,  conventional  line  printers, 
plotting  machines,  and  film  recorders,  Perhaps  the  most  life-like 
pictures  are  continuous  tone  images  produced  on  raster  scan  cathode 
ray  tubes. 

Continuous  tone  display  requires  the  capability  of  defining 
the  light  intensity  of  each  pixel  of  a scan  line  - TV  style.  There 
are  typically  512  scan  lines  per  picture  with  512  pixels  per  line, 
and  256  levels  of  light  intensity  for  each  pixel.  For  a color 
image,  each  pixel  must  know  the  light  intensity  for  each  of  the 
three  primary  colors,  Given  the  intensity  information,  a picture 
can  be  'painted'  pixel  by  pixel,  scan  line  by  scan  line. 

Whereas  the  display  itself  is  strictly  a hardware  problem, 
the  software  problem  is  chiefly  this:  What  intensity  should  each 
of  the  250,000  odd  pixels  have  in  order  to  create  the  desired 
picture?  The  preceding  question  assumes  a microscopic  perspective, 
whereas  the  actual  software  development  proceeds  at  a macroscopic 
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level.  The  overall  software  problem  divides  itself  into  several 

major  sub-problems,  such  as  spatial  orientation  (translation  and 

rotation  of  the  object),  perspective,  hidden  surface  removal,  and 

reflectivity.  This  brief  overview  omits  discussion  of  the  solution 

to  these  problems,  but  the  reader  is  referred  to  a sampling  of 

12  3 

literature  addressing  these  problems.  ’ 

One  point  must  be  made  here,  however.  Continuous  tone 
graphics  concerns  itself  with  surfaces  - specifically  surfaces  of 
mathematical  models.  Consequently,  only  surface  definitions,  as 
opposed  to  line  or  point  definitions,  can  be  used  as  input  data. 

One  way  to  define  an  arbitrary  surface  is  to  approximate  it  as  a 
network  of  discrete  polygonal  elements  (triangles  and  quadrilaterals) 
which  are  defined  first  by  vertices  in  3“D  space,  and  further  by  a 
connecting  perimeter.  Such  a definition  will  herafter  be  referred 
as  a panel  definition. 

The  continuous  tone  pictures  in  this  thesis  were  photographed 
off  a Comtal  Image  Generator,  The  display  files  were  generated  using 
MOVIE. BYU  - a powerful  graphics  package  written  by  Dr.  Christiansen 


1 

* Henry  (1.  Christiansen,  "Applications  of  Continuous  Tone 
Computer-Generated  Images  in  Structural  Mechanics,"  Structural 
Mechanics  Computer  Programs  - Surveys,  Assessments,  and  Availability, 
University  Press  of  Virginia,  Charlottesville,  Virginia,  June  197*1, 
pp.  1003-1015. 

2Henry  N.  Christiansen,  "MOVIE. BYU  - A General  Purpose 
Computer  Graphics  Display  System."  Proceedings  of  the  Symposium  on 
Applications  of  Computer  Methods  in  Engineering,  University  of 
t Southern  California,  Los  Angeles,  August  1 y / / - 

William  M.  Newman  and  Robert  F.  Sproul l , Principles  of 
Interactive  Computer  Graphics  (New  York:  McGraw-Hill,  1973) 
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(of  Brigham  Young  University)  and  Dr.  Stephenson  (now  at  the  Univer- 
sity of  Arizona).  This  thesis  focuses  on  generating  panel  definitions 
from  contour  data  in  a format  compatible  with  the  requirements  of 


MOVIE. BYU. 
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Chapter  3 

CONVERTING  CONTOURS  INTO  PANELS 

A LIMITED  TRIANGULATION  ALGORITHM 

A contour  line  can  be  viewed  mathematically  as  the  inter- 
section of  an  arbitrary  surface  and  a plane.  In  topography,  the 
plane  is  generally  horizontal  at  a specified  elevation.  If  the 
surface  is  closed,  its  contour  lines  will  likewise  be  closed  loops. 

A set  of  contour  lines  on  evenly  spaced  parallel  planes  comprise 
a contour  definition  of  a surface. 

Contour  lines  of  an  irregular  surface,  such  as  found  in 
nature,  do  not  lend  themselves  to  curve  fitting,  or  other  attempts 
at  precise  mathematical  description.  The  most  convenient  numerical 
description  of  a contour  line  is  perhaps  one  where  the  line  is 
approximated  as  a string  of  straight  line  jegments.  This  digitized 
contour  line  offers  two  pieces  of  information:  nodal  coordinates,  and 
connectivity  of  nodes.  Connectivity  is  implied  by  the  sequence  in 
which  the  nodes  are  listed. 

Triangulation  - the  process  whereby  a panel  definition  of 
triangular  panels  is  extracted  from  a contour  definition  - is 
greatly  facilitated  by  observing  the  connectivity  inherent  in  cont- 
our data.  That  connectivity  leads  us  to  explicitly  note  an  obvious 
rule  in  triangulation:  If  two  nodes  of  the  same  contour  are  to  be 
defined  as  nodes  of  the  same  triangle,  they  must  neighbor  each 
other  on  their  contour  line. 
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Also,  no  more  than  two  vertices  of  any  triangle  may  be 
recruited  from  the  same  contour  line  (except,  of  course,  in  the  spe 
cial  case  where  the  entire  area  enclosed  by  that  contour  is  to  be 
capped  off) . 

Triangulation  is  most  logically  carried  on  between  pairs  of 
adjacent  contour  lines,  Consider  this  pair  of  contour  loops  T 
(top)  and  B (bottom). 


Figure  1 

Contour  Pair  Prior  To  Triangulation 

Two  requirements  must  be  met  before  triangulation  commences 
First,  both  loops  must  run  in  the  same  rotational  direction,  and 
second,  the  first  nodes  of  each  loop  must  be  proximate.  Both  rules 
are  met  by  these  loops,  and  they  are  ready  for  triangulation. 

Perhaps  at  this  point  discussion  might  best  center  on  the 
finished  product. 


Figure  2 

Triangulated  Contour  Pair 

Observe  in  figure  2 the  triangulated  contour  pair.  If  one 
were  to  ask  oneself  "How  could  a computer  algorithm  be  taught  to  do 
this?",  a few  ideas  would  assert  themselves.  First,  each  contour 
segment  can  be  considered  to  be  the  base  of  a triangle,  with  the  third 
vertex  being  a node  from  the  other  contour.  Secondly,  each  triangle 
appears  to  be  as  fat  as  possible.  That  is,  the  third  vertex  is 
always  very  near  its  counterparts  on  the  other  contour  line. 

With  these  ideas  in  mind,  consider  again  the  untriangulated 
loops.  Referring  to  figure  3>  triangulation  commences  by  defining 
diagonal  It-lb.  Since  contour  connectivity  requires  1 1 -2 1 and  1b-2b 
as  bases  of  triangles,  there  are  exactly  two  candidates  for  the  first 
triangle:  1t-lb-2t,  and  1t~lb-2b.  Glancing  back  at  the  solution, 


r 
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Figure  3 

Commencing  Triangulation 

it  is  seen  that  1 1 - 1 b-2b  was  selectedv  Moving  on,  once  again 
there  are  exactly  two  possibilities  for  the  second  triangle: 
lt-2b-2t,  and  lt-2b~3b.  This  time,  triangle  lt-2b-2t  is  selected. 
Notice  that  in  each  case,  there  are  only  two  triangles  to  decide 
between,  and  that  the  triangle  with  the  shortest  diagonal  is  chosen. 
This  procedure  continues  until  both  loops  have  been  traversed. 

This  "shortest  diagonal"  algorithm  is  very  easily  implemen- 
ted, and  works  fine  as  long  as  the  two  loops  are  mutually  centered 
and  are  of  reasonably  similar  size  and  shape. 

MAPPING 

The  basic  "shortest  diagonal"  algorithm  fails  for  mildly 
complex  cases.  A typical  example  is  found  in  this  pair  of  offset 
contours. 


Figure  A 
Fa i lure  Example 

Here,  the  shortest  diagonal  search  results  in  a cone. 
Rather  than  abandoning  the  algorithm,  let's  consider  modifying 
the  contour  loops  to  make  them  more  acceptable.  As  mentioned, 
the  algorithm  prefers  contour  pairs  to  be  mutually  centered,  of 
similar  size,  and  of  similar  shape.  The  first  two  requirements 
can  be  met  by  mapping  the  loops  onto  a unit  square  prior  to 
triangulation.  (Mapping  also  tsnds  to  make  the  shapes  more  uni- 
form, though  not  always  enough.  This  problem  is  addressed  in 
the  next  section.) 


Figure  5 

Window  Parameters 


2.  Calculate  A X, A Y,  X,  and  Y. 

3.  Map  onto  a unit  square  centered  at  (0,0)  by  translating 
and  scaling  the  contour  such  that  its  window  matches  the  unit  square's 


window.  The  equations  for  this  are: 

x7=(x-x)/ax 

y'=(y-y)/ay 


Figure  6 

Mapped  Contour  Pair 

With  both  contours  thus  mapped,  they  are  easily  handled  by 
the  original  algorithm. 

A fringe  benefit  of  mapping  is  that  the  resulting  triangles 
tend  to  align  themselves  with  diagonals  that  are  biased  in  the 
direction  of  the  offset.  This  creates  a desirable  longitudinal 
texture. 


ULTIMATE  AMBIGUITY 

A set  of  contour  lines  contains  the  following  mathematical 


information: 


1 
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1.  Exact  coordinates  of  some  points  on  the  surface. 

2.  Approximate  gradients  in  the  X-Y  contour  plane. 

3.  A general  idea  of  the  range  of  possible  Z-gradients, 

A panel  definition  contains  items  1 and  2 and  improves  on 

item  3 by  pinning  down  approximate  Z components  of  surface  gradients. 
Consequently,  there  is  a degree  of  ambiguity  inherent  in  the  tri- 
angulation problem. 

When  two  loops  are  similarly  shaped,  the  ambiguity  is 
negligible.  To  illustrate,  consider  these  two  solutions  of  the 
same  triangulation  problem; 


Figure  7 

Synonimous  Tr iangulat ion  Interpretations 

Since  these  two  solutions  are  different,  one  of  them  is 
probably  a more  exact  approximation  of  the  actual  surface.  But, 
since  the  true  surface  gradients  are  not  available  for  comparison, 
and  since  the  two  solutions  are  so  similar,  either  solution  is 
probably  adequate.  After  all,  contour  lines  form  a skeletal 
framework  that  cast  rather  rigidly  the  shape  of  the  surface. 
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However,  as  the  respective  shapes  of  a contour  pair  become 
increasingly  divergent,  the  ambiguity  becomes  increasingly  pro- 
ncuiicad.  The  following  convolution  provides  a good  example: 


Figure  9 

Non-Synon imous  Tr iangulat ion  Interpretations 

Here,  the  variation  in  interpretation  is  not  as  tolerable. 
Both  solutions  are  reasonable,  yet  one  is  wrong.  Clearly,  more 
information  is  required  to  resolve  this  problem. 

There  are  two  ways  to  provide  the  needed  information.  Firs 
one  could  require  the  contour  planes  to  be  close  enough  together 
that  there  is  minimal  variation  between  adjacent  contour  lines. 

This  approach  has  the  advantage  of  tending  towards  an  exact  descrip 
tion,  ar.d  the  disadvantage  of  being  uneconomical. 

The  second  approach  (adopted  in  this  thesis  by  default)  is 
to  request  user  interaction  to  guide  the  triangulation  over  cases 
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of  excessive  ambiguity.  This  is  a more  general  solution  to  the 
problem.  Here,  the  user  is  called  upon  to  resolve  the  ambiguity  with 
his  knowledge  of  the  true  shape  of  the  surface.  For  the  mechanics 
of  how  this  is  implemented  in  the  computer  program,  refer  to  the 
user  documentation  and  the  example  problem. 


BRANCHING 


An  important  feature  of  this  algorithm  is  the  capability  to 
handle  branching.  Consider  this  simple  case  where  one  contour  loop 
branches  into  two: 


Figure  9 

Simple  Case  of  Branching 


One  way  to  handle  this  is  to  respectively  treat  each  contour 
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as  if  it  were  alone,  neglecting  the  other  branch.  The  resulting 
triangulation  would  appear  like  this: 


Drawbacks  are  that  it  is  uneconomical,  and  it  is  unaccept- 


able in  even  mildly  complex  branching  situations, 

A more  economical,  and  more  general,  approach  to  branching 
is  outlined  in  this  thesis.  The  idea  is  to  treat  all  branches  as 
one  continuous  closed  loop  by  introducing  a new  node  midway  between 
the  closest  nodes  on  the  branches  and  renumbering  the  nodes  of  the 
branches  and  the  new  node(s)  such  that  they  can  be  considered  as 
being  one  loop.  The  Z coordinate  of  the  new  node  is  the  average  of 
the  Z coordinate  of  the  two  levels  involved. 

Plan 


E 1 eva  t i on 


As  seen  from  figure  12,  the  new  node  and  Its  immediate 


neighbors  are  numbered  twice  to  give  the  effect  of  one  continuous 
loop.  Triangulation  can  now  proceed  as  normal.  The  scheme  is 
easily  expanded  to  handle  more  than  one  branch. 

Often,  there  are  several  contour  loops  on  adjacent  planes, 
posing  the  problem  of  loop  connectivity.  Which  loops  should  be 
triangulated  one-on-one,  and  which  are  cases  of  branching? 
Judgment,  in  clear  cut  cases,  can  be  made  on  the  basis  of  window 
overlap. 


Figure  13 

Typical  Problem  In  Connectivity 

Here,  T^  and  clearly  go  together,  and  B2  clearly  branches 

into  T2  and  T^.  Window  overlap  is  best  found  by  default:  IF  they 
don't  not  overlap,  they  overlap. 


Figure  1*1 
Overlap  Test 

The  rectangular  windows  definitely  do  not  overlap  if: 


TY  < 

BY  . 

or 

max 

min 

TY  . > 

BY 

or 

min 

max 

TX  < 

BX  . 

or 

max 

min 

TX  . > BX 
m i n max 

On  the  other  hand,  if  all  four  inequalities  are  false,  the 
windows  necessarily  overlap. 

The  algorithm  works  well  for  mildly  complex  cases,  with 
optional  user  interaction  capabilities  to  handle  complex  branch- 
es 
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Chapter  4 

BRAIN  CONTOUR  DATA 

ORIGIN  AND  DESCRIPTION 

The  brain  data  to  be  used  for  example  purposes  in  this 
thesis  has  an  interesting  history.  In  1 9&7 , the  first  of  several 
movies  was  made  of  a human  brain  at  the  University  of  California 
at  San  Diego.  Using  the  process  of  cinemorphology,  an  entire  human 
brain  was  placed  in  a microtome  capable  of  shaving  off  a slice  25 
microns  thick.  After  each  slice,  a frame  of  movie  film  was  shot. 
The  entire  brain  was  sliced  through,  with  each  successive  newly 
exposed  surface  recorded  on  film.  Every  nth  frame  of  the  movie 
was  exploded  photographically  and  outlines  traced  of  each  distinct 
brain  structure.  Figure  15  shows  a cortex  contour.  In  all,  22 
separate  structures  were  recorded.  The  contour  outlines  were  laid 
on  an  acoustic  tablet  and  a graduate  student  (of  course.’)  selected 
appropriate  nodes  with  the  acoustic  pen.  The  nodes  were  then 
digitized  and  recorded. 

This  digitized  data  base  was  accessed  by  a line  drawing 
graphics  package  which  can  produce  real  time  line  drawing  movies 
(in  color!)  on  an  Evans  and  Sutherland  Picture  System, 

Each  contour  plane  is  referred  to  as  a ’'page",  and  there 
are  98  pages  total,  ranging  from  page  3 at  the  top  of  the  cortex 
to  page  100  at  the  bottom  of  the  brain  stem.  The  data  base  is 
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Figure  15 


Contour  Line  of  Brain  Cortex 

massive  - totaling  78,651  nodes.  Table  1 shows  the  number  of  nodes 
per  structure,  as  well  as  their  page  limits.  Pages  are  spaced 
approximately  1/25"  apart,  corresponding  to  a brain  that  is  roughly 
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Table  1 


INDEX  TO  BRAIN  DATA 
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STRUCTURE 

PAGES 

NODES 

1 

Cortex 

3-78 

52,870 

2 

Caudate 

30-61 

1,922 

3 

Ventricles 

30-84 

4,707 

4 

Fornix 

35-57 

1 ,081 

5 

Putamen 

37-54 

3,075 

6 

Thalamus 

38-58 

1 ,248 

7 

Corpus  Callosum 

41-46 

35 

8 

Globus  Pal  1 idus 

43-52 

725 

9 

Hippocampus 

47-66 

1,576 

30 

Hypothal amus 

50-61 

400 

31 

Pineal  Body 

51-55 

92 

32 

Subthalamic  Nucleus 

50-56 

142 

33 

Red  Nucleus 

52-60 

23S 

34 

Brain  Stem 

54-100 

1,960 

35 

Amygdala 

55-63 

395 

16 

Substantia  Nigra 

56-62 

243 

37 

Cerebel lnm 

59-99 

6,800 

38 

Optic  Chiasm 

60-62 

78 

39 

Mammillary  Bodies 

57-59 

87 

20 

Mesopal 1 ium 

3 9-69 

2,385 

21 

Manual 1 lothalamic  Tract 

43-56 

303 

22 

Septum 

42-49 

289 

■! 


A 


23 


( BRAIN  DATA  FORMAT  TRANSMUTATION 

The  brain  contour  data  arrived  at  BYU  on  magnetic  tape  as 
l6bit  integers  in  binary  format.  The  data  are  grouped  into  22 
structures,  which  in  turn  are  divided  into  segments.  A segment 
is  a string  of  contour  points  and  a contour  line  is  formed  from 
one  or  several  segments.  Segments  represent  portions  of  surfaces 
which  are  shared  by  two  structures.  Hence,  a contour  line  that 
is  composed  of  say  5 segments  is  bordered  by  5 neighbors. 

Segment  definition,  which  was  initially  imposed  on  the  data 
to  facilitate  line  drawing  display,  somewhat  hampers  triangulation 
because  all  contours  must  be  reconstructed  from  their  constituent 
segments  before  triangulation  can  commence.  The  problem  is  aggra- 
vated because  the  segments  are  randomly  sequenced  and,  furthermore, 
no  convention  is  observed  in  clockwise  and  counterclockwise  ordering 
of  nodes. 

Segment  definition  is  illustrated  by  this  typical  configura- 
tion where  3 closed  contour  loops  are  defined  by  6 segments: 


2h 

Notation  wise,  Pnj  is  the  first  node  of  segment  n and  Pn2  the  last 
node  of  segment  n.  The  transmutation  algorithm  begins  by  assigning 
segment  1 to  loop  1.  A search  is  made  for  the  nearest  neighbor  of 
Pj2  which  is  P^j,  and  segment  5 is  appended  to  segment  1.  Next, 
the  nearest  neighbor  of  P^2  "> s sought.  Its  nearest  neighbor 
is  P32.  This  indicates  that  segment  3 is  sequenced  in  an  order 
contrary  to  that  of  segments  1 and  5.  Consequently,  segment  3 is 
joined  to  loop  1 in  reverse  order.  Since  P^j  neighbors  P^j,  the 
loop  is  complete.  A flag  is  set  for  loops  1,  3 and  5 preventing 
future  assignment.  This  logic  repeats  until  all  segments  are  joined 
to  a loop. 

i 

It  is  important  to  impose  the  convention  that  loops  run 
uniformly  in  a clockwise  (or  counter-clockwise)  direction.  To 
enforce  this  convention,  all  nodal  angles  of  a contour  line  are 
summed  for  monitoring  In  the  following  manner: 


Figure  17 

Determining  Direction  of  Rotation 
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Theta  is  the  angle  by  which  each  succedent  vector  deviates  from  a 
straight  line.  The  sum  of  all  such  angles  will  be  360  degrees  for 
counterclockwise  sequencing  and  -360  degrees  for  clockwise.  Theta 
is  computed  from  vector  cross  and  dot  products. 


(sin9)  cos9>0. 

9 

(sin0)+9O  , cos6<0  and  s?n0>O. 

(sin0)*9O  , cos0<O  and  sin0<O. 

In  a few  instances,  the  brain  data  invalidates  this  approach  by  having 
a contour  line  cross  itself  like  this: 


sin0= 


M 

AB 


COS0  = 


t-t 

AB 


sin'1 
0=  sin  1 
sin"1 


This  error  causes  a figure  8 which  results  in  10  approaching  0.  Often, 
this  causes  a violation  of  the  rotation  convention. 


Chapter  5 


ECONOMIZING 

NODE  ELIMINATION 

If  a data  base  is  too  refined  (i.e.  contains  nodes  you  could 
do  without)  it  is  desirable  for  reasons  of  economy  to  eliminate  the 
less  essential  nodes.  Consider  this  node: 


The  node  is  accepted  (or  rejected)  upon  the  following 
criteria: 


reject  accept  reject 

node  node  node 


Figure  20 

Node  Reduction  Flow  Chart 
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S . , S and  6 . are  user  definable  parameters.  Every  node  is 
min  max  min 

screened  using  this  logic.  To  assure  acceptance  of  every  node,  all 
three  parameters  may  be  set  to  zero. 

This  is, a logical  place  to  interject  a few  thoughts  on 
interpolation  of  new  nodes.  Since  a digitized  contour  line  is  an 
approximation  comprised  of  a series  of  straight  line  segments,  it  is 
reasonable  to  assume,  and  important  to  prescribe,  that  nodes  are 
selected  such  that  the  digitized  approximation  does  not  deviate 
intolerably  from  the  actual  contour  line.  This  would  imply  a 
correlation  between  nodal  density  and  contour  line  curvature.  That 
correlation  suggests  that  it  is  desirable  to  re-distribute  nodes 
around  the  contour  loop  according  to  a curvature  vs.  node  density 
function  using  curve  fitting  procedures.  This  is  an  appealing  thought, 
since  it  would  reduce  angularity  in  the  continuous  tone  display. 

This  would  be  great,  provided  the  actual  surface  isn’t  angular.  Of 
course,  however,  that  assumption  is  not  always  valid.  Take,  for 
example,  a simple  four  node  contour  definition  of  a square.  If 
curve  fitting  were  imposed  in  an  attempt  to  extrapolate  extra  nodes, 
the  result  would  tend  towards  a circle.  Angularity  is  reduced  at  the 
expense  of  accuracy.  In  conclusion,  the  burden  of  providing  acceptable 
data  rests  with  the  person  who  actuallv  does  the  digitizing. 

QUADRILATERAL  FORMATION 

.'or  economy  of  storage,  it  is  desirable  to  join  pairs  of 
adjacent  triangles  together  into  quadrilaterals.  Not  all  ajoining 
triangles  are  thus  combined.  The  decision  is  made  on  the  basis  of 
how  Tped  the  resulting  quadrilateral  would  be,  that  is,  the  angle 
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by  which  the  two  triangles  are  out  of  plane.  That  angle  is  easily 
found  using  vector  algebra. 


s i n a = 1 - 2 

L1  2 


If  sina<sin  a , the  two  triangles  are  redefined  as  a 
max 

quadrilateral.  a is  user  definable  and  defaults  to  ^5  degrees, 
max 
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Chapter  6 

FORTRAN  IMPL IMENTAT! ON 

A few  explanatory  remarks  are  offered  here  to  the  reader  who 
is  bent  on  deciphering  the  source  code. 

As  mentioned,  the  brain  data  arrived  at  Brigham  Young 
University  in  the  form  of  16  bit  integers  on  tape  in  binary  format. 

To  facilitate  use  on  the  DEC-10  computer,  these  data  were  re-formatted 
into  7 bit  ASCII  data  files,  one  file  for  each  of  the  22  structures, 
with  8 integers  per  line.  The  first  two  integers  of  a file  comprise 
the  "structure  heading",  the  first  integer  being  the  structure  number, 
and  the  second  being  the  number  of  segments  in  the  structure.  These 
two  integers  are  ignored  by  the  tr iangu lat ion  program. 

The  next  four  integers  form  the  segment  header  of  the  first 
segment.  The  first  integer  of  the  segment  header  is  the  page  number 
or  horizontal  level  of  the  segment.  The  page  numbers  range  from 
”3  to  -99.  The  Z coordinate  of  the  segment  is  computed  from  the 
formula  RZ=- (Z+51 )*A50*SCALE  where  RZ  is  the  Z coordinate,  Z is  the 
page  number,  and  SCALE  is  the  scale  factor.  The  second  and  third 
integers  in  the  segment  header  are  ignored.  The  fourth  integer,  NPL, 
is  the  number  of  nodes  in  the  segment. 

Immediately  following  the  segment  header  are  the  X-Y 
coordinates  of  the  segment  nodes,  totalling  NPL  pairs.  The  next 
segment  header  immediately  follows  the  last  node  of  the  preceding 
segment,  so  there  is  an  uninterrupted  string  of  integers  from  start  to 
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finish  in  the  data  file.  Nodes  of  all  segments  of  the  same  page  are 


stored  in  array  P2.  The  segments  are  joined  together  to  create 
closed  loops  (as  described  in  chapter  k)  and  stored  in  array  P3. 
Node  elimination  is  imposed,  and  the  nodes  are  finally  stored  for 
tr i angul at  ion  in  array  P. 

Two  pointers  are  used  in  accessing  the  nodes  in  array  P. 
The  Loop  Pointer  - LPP  - indexes  the  global  loop  numbers  of  the 
first  loop  on  any  contour  level.  Pointer  PI  addresses  the  global 
node  number  of  the  first  node  on  a loop.  For  example,  the  X 
coordinate  of  the  n*"^  node  of  the  loop  on  the  k**1  contour 
level  would  be: 

P((Pl(LPP(k)+(j-l))+n-l)),1) 

Knowing  these  conventions,  the  fortran  coding  should 
be  relatively  lucid.  The  coding,  written  for  use  on  a DEC-ID 
system,  is  given  in.  the  appendix.  An  additional  assembly  language 
Tektronix  interface,  GRTEK.REL,  must  also  be  loaded. 
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DRAIN  DATA?  V 

FILENAME  OF  INPUT  HAT  A?  T :p 
REAIi>  HELP 

PARAMETERS?  TOTALS?  LEUEL?  SCALE?  EXIT 
BRANCH?  MANUAL?  CLOSE?  DEUICE?  KLOCKWISE 


REAIi>  P 


MINIMUM  SEGMENT  RNGLE=  15 

MIN.  & MAX.  SEGMENT  LENGTHS-  .3? a 

READ>  L. 

Z—SPACING=  1 


LEUEL  RANGE-  1 ? CO 

HATH  ENDED  AFTER  LEUEL  3£ 
BRANCH > HELP 


AUTOMATIC?  HARP?  MANUAL?  INSPECT ? 
BRANCH > T 


SINGLE?- CAP?  EXIT  ? TOTALS 


955  NODES.  1 ELEMENTS 
BRANCH > A 


START  WITH  WHICH  LEUEL?_1 
POST-EDIT?  V 


Figure  23 
Computer  Dialogue 

selected.  We  opt  to  set  the  minimum  segment  angle  to  15  degrees,  S . 

m i n 

to  .3  and  S to  2.  (These  parameters  are  explained  on  page  26). 
max 

Having  returned  to  the  READ>  prompt,  we  now  choose  to  read  in  the 
data,  being  satisfied  with  the  default  values  for  SCALE,  CLOSE,  and 
KLOCKWISE.  After  setting  the  LEVEL  parameters,  the  computer  goes  to 
work  reading  in  contour  segments,  reconstruct ing  them  into  loops. 


thinning  them  out,  and  assuring  that  all  loops  run  in  a clockwise 
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direction.  The  algorithm  encounters  end  of  file  before  60  levels  are 
read  in,  and  informs  us  that  all  32  available  levels  have  been  read  in 
and  processed.  The  TOTALS  command  receives  the  response  that  there 
are  955  nodes  in  array  P,  which  means  that  over  half  of  the  available 
1922  nodes  have  been  thinned  out  by  the  node  elimination  algorithm. 

The  Caudates  are  quite  regular,  requiring  little  or  no 
interaction  to  triangulate  properly,  so  the  AUTOMATIC  option  is 
invoked  beginning,  as  usual,  with  level  1.  As  a safeguard,  post- 
editing is  requested.  If  the  data  were  unquestionably  obedient,  the 
post-editing  could  justifiably  be  circumvented,  but  this  way  good 
results  are  guaranteed. 

The  algorithm  now  proceeds  to  first  determine  connectivity, 
then  to  triangulate  all  loops  implicated  in  the  the  window  overlap 
connectivity  check,  and  finally  to  display  the  resulting  panel 
definition.  The  user  glances  at  each  successive  display  and  grants 
acceptance  with  a carriage  return,  or  occass ional ly  rejects  the 
triangulation,  as  the  case  may  be.  (Usually  it  is  immediately  clear 
when  triangulation  is  unacceptable.  Normally,  failure  occurs  in 
areas  where  the  two  loops  are  excessively  dissimilar,  and  the  resulting 
panels  are  often  bizzare). 

Two  of  several  simple  one-on-one  triangulations  are  shown 
in  figures  25  and  26  , each  of  which  is  accepted.  However,  the 
branching  loops  in  figure  27  are  triangulated  incorrectly,  and  a 
change  is  requested.  This  change  is  granted  through  erasure  of  the 
screen,  re-drawing  of  the  untriangulated  loops,  and  issuing  of  the 
TRIANGULATE>  prompt.  The  INTERACTIVE  command  is  given.  As  always. 


only  the  first  letter  of  the  command  is  required.  Referring  to 
figure  28  the  nodes  are  numbered  for  identification  - even  for  the 
top  loop  and  odd  for  the  bottom. 

The  INTERACTIVE  command  allows  the  t r i angu la t ion  to  be 
controlled  by  allowing  specification  of  nodal  delimiters  between 
which  t r iangulat ion  will  occur.  For  example,  if  the  delimiteres  1,1 
for  the  top  and  1,2  for  the  bottom  were  chosen,  only  one  triangle 
would  be  formed.  Basically  the  selection  of  delimiters  is  a trial 
and  error  process.  The  user  delimits  as  large  a span  as  reasonable. 
If  the  resulting  t r i angu 1 a t i on  is  adequate,  great!  If  not,  try  again 
with  a smaller  span.  A general  rule  might  be  to  procede  quickly  over 
sections  where  the  two  contours  are  similar,  and  cautiously  over 
sections  where  they  are  dissimilar. 
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Referring  to  figure  26,  the  triangulation  was  successful 
through  node  17  of  the  bottom  loop,  but  then  failed  to  traverse  the 
branch  properly.  This  suggests  the  selection  of  delimiters  1,28 
for  the  top  and  1,17  for  the  bottom.  The  resulting  tr i angul at  ion 
shown  in  figure  28  looks  good.  Now,  one  more  span  (top:28,33  ; and 
bottom: 1 7,20)  should  suffice,  and  figure  29  confirms  the  hope.  Now 
that  the  entire  circuit  is  complete,  AUTOMATIC  mode  is  re-entered  and 
triangulation  proceeds  smoothly  to  the  conclusion.  Upon  exiting  from 
the  program,  the  panel  definition  is  written  into  a disk  file, 
available  for  display  using  MOVIE. BYU.  A continuous  tone  image 
created  by  that  panel  definition  is  shown  in  figure  24. 
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Figure  29 
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Chapter  8 

PICTURES 

This  chapter  presents  examples  of  the  finished  product.  It 
is  difficult  to  judge  how  true  to  life  the  images  are,  due  to  the 
highly  esoteric  nature  of  the  subject  matter.  Nonetheless,  it  is 
generally  evident  that  the  triangulation  algorithm  has  performed 
reasonably. 

The  first  structure  presented  is  the  Brain  Stem.  The  Brain 

Stem  was  a straighforward  triangulation  problem.  There  are  no 

branches,  and  there  is  no  sreious  variation  in  the  shape  of  its 

respective  contours.  The  only  difficulty  was  a case  of  illeagle 

rotational  direction,  as  described  on  page  25.  Other  than  that, 

the  entire  triangulation  was  handled  automatically.  Two  panel 

definitions  were  generated  for  the  brain  stem.  The  first  file 

was  made  wi th  S . =S  =0  . =0,  Z-SPACING=1,  and  LEVEL  RANGE  =5,36. 
mm  max  min 

The  resulting  panel  definition  had  1609  nodes  and  1991  panels. 

The  second  brain  stem  panel  definition  has  378  nodes  and  1*92  panels. 

It  was  generated  wi th  S . =.2,  S =1,0.  =15.  Z-SPACING  =k,  and 

min  max  min 

LEVEL  RANGE  =5,1*0.  Line  drawings  are  shown  in  figures  30  and  31, 
and  continuous  tone  images  in  figures  3 4 and  35. 

The  next  pair  of  images,  figures  36  and  37  are  detailed 
studies  of  the  thalamus.  Here,  1*  loops  branch  into  2,  then  2 into  1. 
Figure  37  has  a more  biological  look  due  to  Gouraud  smooth  surface 
simulat ion . 


k2 


Figure  38  is  a strik  i ng  composition  of  6 different  struc- 
tures, each  in  proper  relative  orientation.  The  sructures  are  identi- 
fied in  figure  32  and  figure  33  shows  a line  drawing.  To  enable 
so  many  parts,  the  larger  structures  have  coarser  panel  definitions. 

All  of  the  preceding  examples  were  triangulated  with  virtually 
no  user  interaction.  The  cortex  slice,  in  figure  39  was  not  so 
oblidging.  This  image  is  presented  to  demonstrate  the  degree  of 
complexity  the  algorithm  can  accomodate.  To  help  orient  the  reader, 
this  image  represents  about  a £"  thick  slice  of  the  cortex,  centered 
about  i"  from  the  top  of  the  brain.  The  many  oddly  shaped  holes  are 
due  to  the  fact  that  the  top  most  convolutions  are  decapitated  in 
this  view.  This  panel  definition  - shown  here  with  smooth  shading  - 
consists  of  1752  nodes  with  1778  panels.  This  particular  data  did 
not  cooperate  with  the  automatic  algorithm,  and  required  nearly  2 
hours  time  to  interactively  triangulate.  The  difficulty  was  not 
so  much  the  complexity  of  the  shapes,  but  the  dissimilarity  between 
adjacent  contours.  Also,  the  data  had  a disproportionate  number  of 
glitches.  Nonetheless,  the  continuous  tone  image  is  quite  convincing. 

One  final  image  is  the  panoramic  shot  of  Mt.  Timpanogos  as 
it  might  be  seen  from  an  airplane  flying  to  the  west  of  Timp. 

This  is  submitted  to  illustrate  a possible  application  of  triangulation 
to  display  of  topographic  surfaces. 
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Figure  35 

Brain  Stem  with 
*♦92  Panels 


9 Figure  36  Figure  37 

Thalamus  with  Thalamus  with 

Flat  Shading  Smooth  Shading 


Figure  3** 

Brain  Stem  with 
1991  Panels 
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CONCLUSIONS 

This  thesis  purports  to  present  a general  solution  to  the 

problem  of  converting  a contour  definition  of  an  arbitrary  surface 

into  a panel  definition.  That  assertion  is  rigorously  tested  by 

the  brain  data,  and  experience  with  that  highly  complex  data  base 

lends  credence  to  the  claim  of  a general  solution,  Total  user 

interaction  capabilities  virtually  gaurantee  a general  algorithm. 

Work  might  be  done  on  reducing  the  amount  of  user  dependence 

in  the  algorithm,  though  most  reasonable  cases  require  no  interaction 

at  all.  Also,  it  would  be  helpful  to  improve  graphical  interaction 

by  using,  for  example,  a tablet  to  input  interaction  parameters. 

Study  might  also  be  made  on  how  the  economy  parameters  (S  . , S , 

m'n  max 

0 . , and  a ) effect  the  continuous  tone  image, 
min  max 


i 
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APPENDIX  A 


COnPUTEk  PKOliKAM 


CO«HCN/C/NL (?).LPSTk (?,51.NPT,nJ,0Z,BD,TEKT.PE. 

1 FX(100.4),P1  («OO),TRir,IPL(2,5,10).NIPL(S),lPLI»MP*,*K.lP(a#?000)' 

COP«CR/KOOF/P?(nOO,3).P3tllOO,3) 

COmmcn/7EK/XMAC,Y«AW.SF,P(?000,3),IY02 

INTEGER  LPP(«0).Pl,PIPSI,nPI<S).t«PlC(5).KPL*f2.5) 

INTEGER  LT(<IO),C(6,e).NPL(flO).NPLN(40),NCONC(a0),SURR{B,O0) 

REAL  E(«).ID(8) 

LOGICAL  FLAG(oC).CLOSE,DATAF,TEkT,BD,CR,PE 

INTFGFR  ?S,PP.SP,EP.DP»CR.fiC,TO,BO,>>PANCH{JO),STACK(lG; 

OATA  SAS,SmINT,Smaxt,RANG,07.SCalE/.3,.01,1.,.7,.«5,.0001/ 

OAT  A 7S.KP, L. OP. LIC. PM .FP1 . IZP/J ,0, 1 .2, 1 ,0, I ,0/ 

0*TA  NLA.  NLP.  NL  V.  NP.LTCl). PHI), 1Y  02/0, 0.0. 0.0. 1.0/ 

OATA  LPP(1).IPLI,nPm»X/1,1,2/ 

OATA  TE*  T,  CLOSE.  C*TAF,eD,CR/.  TRUE.,.  TRUE.,.  TRUE.,  .TRUE.,  .TRUE.’/ 

OAT*  K‘ I PL /l, 0,0, 0.0/ 

IPLd.I.DsJ 

SmIN=(CZ»SmINT)*«2. 

SmAXsCCZ*5MAXT)*a2. 

C INPUT  specs 
TYPE  1 

1 F OBM A T ( • BRAIN  DATA?  *,») 

ACCEPT  10.ANS 
IF(anS.E0.*Y*)G0  TO  2 
B0=. FALSE. 

SAS=0. 

SMINaO. 

Sw  A X r I 0 
SCALE=1 . 

2 TYPE  3 

3 FORM*? ( 1 FILENAME  OF  INPUT  DATA?  »,S) 

ACCEPT  O.INAME 

4 F0R“AT(A5) 

0PEN(UNITr21,FILEsINAME) 

C INITIALIZE  OISK 

»EA0(2I,5)A,P,TZ1 
BACKSPACE  21 

5 FORMATOG) 

C ACCEPT  Rt  * D COMMANDS 
GO  TO  B 

6 TYPF  7, NLA 

7 FORMAT!*  LAST  LEVEL  HAS  *.I2) 

0 TYPE  9 

9 F0OMATC  READ*  '.*) 

ACCEPT  10, RC 

10  FORMAT(Al) 

IF (PC  .FG  . * B * )GC  TO  76 
lF(RC.EG.'P*)Gr  TO  20 
IFLPC.EC.'L'.OR.RC.EO.'H'JGO  TO  25 
IF(PC.FC.'S')GC  TO  IP 
IF (PC.FO. 'T * 1TYPE  24.NJ.NPT 
IF(RC.EC.  *E ' )GC  TO  139 
IF(PC.EC.*K* )GC  TO  16 
IF(RC.EC.'C*)GO  TO  1« 

IF (RC.EC. *0* >GC  TO  12 
IF(RC.EC.*T*)GC  TO  S 
TYPE  11 

11  F OPma T ( * PARAMETERS. TOTAl S. LEVEL. SCALE, EXIT, BRANCH, 

1 P amial, CLOSE. PE  V I CE.K LOCK  RISE ' ) 

GO  TO  6 
C SET  DEVICE 


50 
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1 


51 


12 

IS 


c 

10 

15 


C 

lb 

17 


C 


16 

19 


C 

20 

21 


22 


C 

20 

C 

25 

26 


27 

C 

26 


TYPE  13 

FORMAT!'  tektrcnix  SCOPE?  ',$) 

ACCEPT  10, ANS 
TEKTs. FALSE. 

1F( ANS.EO.  'Y' )TEKTs. TRUE. 

60  TO  6 

SET  'CLOSE'  FLAG 
TYPE  15 

FORMAT!'  CLOSE  ALL  LOOPS?  *,S) 

ACCEPT  10#  ANS 
CLOSES. FALSE. 

IF (ANS.EO.  • Y' ) CLCSE  = . TRUE. 

SET  CLOCKWISE  FLAG 
GO  TO  6 
TYPE  17 

FORMAT!'  CLOCKWISE  0RDERING7  ',*) 

ACCEPT  10, ANS 
CU=. FALSE. 

IF  ( ANS.EO.  *Y'  Kws.  TRUE. 

GO  TO  16 

SET  SCALE  FACTOR 
TYPE  19 

FORMAT!'  SCALE  FACTOR®  '#*) 

ACCEPT  *, SCALE 
0?sZS*a50.*SCALE 
SMINr (P7*SMINT)**?. 

SMAXsI0Z*SMAXT)**2. 

GO  TO  8 

SET  PARAMETERS 
TYPE  21 

FORMAT ( ' MINIMLM  SEGMENT  ANGLES  *,S) 

ACCEPT  * # S A N G 

SAS=SIND(SANG) 

TYPE  22 

FORMAT!'  MIN.  A MAX.  SEGMENT  LENGTHS:  ',*) 
ACCFPT  *,SMINT,SMAXT 
SmIns!P2«sminT)**2. 

Sma*s!0Z*SmaxT)**2. 

CO  TO  8 

totals 

FORMAT ! tH  ,10,'  NODES  ',10,.'  ELEMENTS') 

SET  SPACING  AND  RANGE 
TYPE  26 

FORMAT!'  Z-SPACINGs  ',$) 

READ!5#*#END=8,ERRs8)ZS 

DZs.05*ZS 

IZPsZS-1 

TYPE  27 

FORMAT!'  LEVEL  RANGES  *,S) 
READ!S,»,ENDs8,ERRsfl)NLS#NLF 
LINDsNLVAl 
INITIALIZE 

READ!2l,5,END=fc9)IID!J), Jsl,DP),TZ 
6ACKSPACF  21 
IF !DP.EC.6)RACKSPACE  21 
IF!TZ1.E0.TZ)GC  TO  29 


NLAsNLA*! 

IFfNLA.GE  .NLP ) GO  TO  31 

TZlsTZ 

IZPsIZPAl 


L 


52 


IF((I7P.GE.ZS).AND.(NLA.GE.NLS))G0  TO  31 
NPr  0 
KPr  0 

29  KPsKP+1 

C READ  FRO”  RISK 

READC21,5,FND  = 69HI0CJ) , J=1  »DP)  ,Z,  A,DUM,NPL(KP) 
I . (((P?(J+NP,I)),lsl,2),JslfNPL(KP)) 

NP=\P+NPL (KP) 

LT(KP  + n=LT(KP)  + NPL(KP) 

BACKSPACE  21 

C SET  TATA  FILE  POINTER 

DPr2»NPL (KP)+DP+12 
SO  OP=DP-B 

IF(PP.GT.6)G0  TO  30 
GO  TO  26 

C CHECK  SEGMENT  INTERCONNECTIVITY 

31  LCs 1 
NP3sP 

IF(RC.NE.'M')GO  TO  37 

C GRAPHICS  FOR  MANUAL  RE-CONSTRUCTION  cf  loops 
C FIND  WINDOW 

E ( 1 )sP2( 1 , 1 ) 

E(?)SEC1) 

E(3)sP2(l#2) 

E(<I)«EC3) 

DO  32  I s2 # NP 

IF(P2(I,1).LT.E(1))E(1)sP2(I,1) 
IF{P2(I,l).GT.FC2))EC2)sP2CI»n 
IF(“?(l.2).LT.t (3))E(3)sP2(I»2) 

32  IF (P2(I*2).GT.F (R))E (P)=P2(I»2) 

DX=E(2)-E(1> 

0Y=EC4)-E(3) 

XPARcCE  C2J+E  Cl))/2. 

YRAPs(E(«)*E(3))/2. 

SF  = OX 

IFfOr.GT.SFJSFeOY 

SF=700./SF 

C PAINT  * LABEL  SEGMENTS 

IF(.NOT.TEKT)GO  TO  37 
CALL  BOX 
DO  36  1 3 1 * KP 
INDlsLTCI)4l 
IN023LT(I+1) 

IXIrSF*(P2(INDl,n-XBAR)463« 
lY1=SF»(P2(IN01*2)-YBAR)+390 
CALL  MVTOflXl.IYl) 

00  33  J=IN01#IN02 

IXrSF*CF2(J,l)-XRAR)+63« 

Iy=SF«(P2CJ,2)-YBAR)4390 

33  CALL  VCTCCIX.IV) 

IX  = SF«(P2(I NO  1 + 1,1 )-XBAR)+62« 

IY=SF+(P2(IND1*1 ,2)-YPAR)+3S0 
CALL  MVTC(IXflV) 

CALL  mVT0(1X,IY) 

CALL  ALMODE 
TYPE  3«* I 

3«  FnRM*T(lH+tJ.I2.,l') 

TXr$F  + (P?(lN'0?-l*l  )-XRAR)+N2« 
IY=SF»(P2CIN02-1,2)-YPAR)+380 
CALL  MVTO(IX.IY) 
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CALL  ALPODE 
TVt>E  35,1 

35  FOBmAT(1H»,S,I2.'2') 

36  CONTI  NL'E 

CALL  MV  7 0 ( 0 , 767 ) 

call  alpode 

37  IF(<P.EC.l.AND.RC.NE.'P')GO  TO  51 

IFC.MOT.POJGO  TO  51 

DO  38  1=1, KP 
36  FLAG(I)r. FALSE. 

DO  50  1L=1.XP 
IF(FLAG(IL))GO  TO  50 
JPC  = 1 
JLC=IL 

DO  39  1=1, NPL ( TL ) 

NP3=NP3*1 
DO  39  J= 1 , 2 

39  p3(np3»J)=p?CLT(IL)+I,J)*SCALE 

Npl  (LC)  = NPL(IL) 

C F 1 PD  Th=  CLOSEST  ENDPOINT 

00  IFtBC.NE.’H’lGC  TO  03 
JP=2 

IF(JPC.E0.2)JP=1 
T YPF  01,JLC,JP 

01  FORMAT!'  JOIN  ',211,'  TO  ',$) 

ACCEPT  02,JLC,JpC 

02  F OR- A T (211) 

IF( JLC.KE.OIGO  TC  05 
GO  TO  50 

03  Dh=iFa35 
*1=P3(NP3, 1 l/SCALE 
Y1=P3(AP3,2)/SCALE 
DO  00  JL  = 1 , KP 

IF (FLAG! JL ) )G0  TO  00 
DO  ou  JP=1,2 

IF((JL.EO.IL).*NO.(JP.E0.2))GO  TO  00 
JPRr JP 

IF(JP.EC.21JPRsNPL(JL) 

DIST=(Xl-P2(LT(JL)*JpH,l ))**2.*(T1-P2(LT{JL)*JPR,2))**2, 

IF (DIST.GT.OM)GO  TO  00 

JLCsJL 

JPC=JP 

OMsOIST 

00  CONTINUE 

05  FLAG(JLC)=.TRUF. 

IF( JLC.NE.IL) GO  TO  06 

LC=LC*1 

GO  TO  50 

06  IF(JPC.E0.2)G0  TC  08 
DO  u7  1=1 ,NPL ( JLC ) 

NP3=NP3+1 

00  07  J= 1 , 2 

07  . P3(NP3. J)=P2(LT(JLC)AI, J)«SCALE 

NPL (LC)=NPL(LC)*NPL(JLC) 

GO  TO  00 

08  DO  09  1 = 1, NPL  (JLC  ) 

1 I=npl ( JLC )♦!•! 

DO  09  J= 1 , 2 

09  P3(nP3aII,J)=P2(LT(JLC)aI,J)*SCALE 
NP3=NP3*NPL(JLC) 
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npl(LC)=NPLCLC)+npI(JLC) 

GO  TO  00 

50  CONTINUE 
GO  TO  53 

51  00  52  1st, NP 
00  52  J=l,2 

52  P3(I,J)sP2(I,J1»SC*LE 
LC=KP+l 

C RESFT  FLAGS 

53  LCsLC-1 
NLVrNLV+1 

LPP(NLV*n=lC4lPP(NLV) 

C THINNING  OUT 
FIBSTsO 

BZ=-(Z*51)*«50*SCALE 
IF(.NOT.PD)RZ=Z*SCALE 
00  67  Hs 1 , LC 
RP=1 

NLPsNLP+1 
SP=2 
EPsl 
OCsO. 

N=P1 (NLP)-l 

P(N+l,1)sP3CFIPST+I,l) 

P(n*1,2)sP3(FIRST+1,2) 

P(N*I,3)=R7 

C INITIALIZE  N I NOON  MATRIX 

EX (NLP, ! )=P3(FIRST  + 1, I) 
f X(NLP,2)=P3fFIRST+l, 1) 

EX(NLP,3)=P3(F1RST*J,2) 

EX(slp.o)=P3(FIRST*1,2) 

IF(npl(m).GT.2)G0  TO  50 
PI (NLP+1 )=P1 (NLP)ANPL(M) 

IF ( npl (M).EG.IIGC  TO  67 
P(N*2,I)=P3(FIRST+2,1) 

P(N*?,2)=P3(FIRST»?,2) 

P(n*2,3)=RZ 

50  XA  = P3(FIR?T  + RP*EP.n-P3tFIRSTARP,l) 

YA=P3(PIPST*op*tR,2)-P3(FlRST+RP,2) 

A|  sXA*XA«VA*VA 
IF (AL.LT.SMIN)GO  TO  56 

55  X«rP3(FIPST*RP4EF*l, 1 I-P3(FIRST,RP4EP,  I) 

VR  = B3tF  IP$T*bp4EPM,2)-P3(F  IRSTaRP  + EP,2) 
BL=XP»XP*YB*VB 
IF(bl*AL.EO.O.)GC  TO  56 
ST  = (XA*YP-XB<*YAl/SORT(AL<»BL> 

TF(ST.LT.-l.)STa-l. 

IFtST.GT.l.ISTsI. 

IF(ABSfST).GT.SAS)GO  TO  57 

C ANGLE  OP  SEGMENT  LENGTH  IS  TOO  SMALL.  ELIMINATE  THE  NODE. 
IF ( ( (XA4X9)»*2,A(YA4YF)«»2.).GT.SMAX)G0  TO  57 

56  EP=Ep»l 
IFf (RP+EP).EO.KPL(H))GC  TO  5P 
GO  TO  50 

C ACCEPT  THE  NODE. 

57  P(N4SP,1)sP3(FIBST4RP*EP»1> 
p(n*SP#2)=P3(FIRST4RP4EP,21 
p(n*SP,3)=RZ 

C SUM  ANGLES  TC  DETERMINE  DIRECTION  OF  ROTATION 
OC*OC»ASIN(ST) 

■ 
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IF((XA«xP+YA*yP).GF.O.)GO  TO  571 
IF(ST.LT.0.)DC*DC-1.5708 
IFCST.GT.O.)DC*OCy1.5708 
C DETERMINE  LOOP  EXTRFPES. 

571  IF(P(N*SP, 1 ) .LT ,FX(NLP. 1))EX(NLP,1 )*P(N»SP,1) 

IF(P(n*sp.  1 ) .GT.EX(VlP,2))EX(KiLP,?)=P(K»SP,  t) 

lF(PfM+«P,?).LT.rx(NLP,3))EX(NlF,3)=P(N*$P,2) 

IF(P(K'  + EP,2).GT.EX(NLP,«))tX(NLP»«)=P(MtSP,2) 

X AsXB 

YA  = YB 

*L*BL 

PPsBP+EP 

SPsSPy 1 


C 

58 


60 

61 


I 


63 


60 

C 


65 


66 

67 
C 

68 


69 

70 

71 
C 


JF(PP.EO.(K'PLCP)-mGO  TO  58 
E P=  1 

GO  TO  55 
ENDPOINTS 

N S P = 5 P 
00  59  1=1,2 

IF(P3(FIRST  + NPt(Y),I).KE.F(N*l,n)G0  TO  60 
NSP=NSP-1 
GO  TO  62 
00  61  1*1,2 

P(N*SP. I )*P3(FIRST*NPl (K), I) 

P(N*SP,3)*RZ 

IF(.MOT.CLOSE)GO  TO  60 

NSPsNSPYl 

00  63  1*1,2 

P(K'  + NSP,  I )=PCN+1, 1) 

P(N»KSP,3)*9Z 
PI (NLP*i )=P1  CMP) »NSP 
INSURE  CLOC**>ISt  ORDERING 

!F(OC.LT.O..»NC.CN)GO  TO  67 

IF ( (DC.GT.O.  ) ,»ND. (,N0T.CW))G0  TO  67 

DO  65  1*1, MSP 

00  65  J* 1 , 2 

p?(I,J)=P(NtI»J) 

00  66  1*1, NSP 
DO  66  J*  > » 2 
II=NSP*1-I 
P(5*II,J)*P2CI,J) 

F IMSTrF IRST+NPL (P) 

RESET  POINTERS 
KP  = 0 
IZPsO 
NP  = 0 

IF(NLA.LT.NLF)GO  TO  28 
GO  TO  71 
Type  70, NLA 

F OR"AT  f 1 data  E NCEP  AFTER  LEVEL  ’,12) 
OATAFs. FALSE. 

N JsN.NSP 


OETERPINE  CONCENTRICITY 
DO  75  p*  1 , NL V 
NCONC (P)*0 
L 1 1 *LPP (P) 
LI2=LPP(PY1 )-l 
DO  75  I*LI1,LI2 
00  75  J*L  I 1 , L 12 
IF(I.E0.J)G0  TC  75 


J 


c 

72 

73 
C 


7<i 


75 
C 

76 

77 


70 

C 

79 

00 


C 

01 

02 


03 


CMFCK  FOP  TRIVIAL  REJECTION 
DO  72  Ks 1 , 3, 2 

IF(EX(J,K*1).LT.EX(X.k))C0  TO  75 
IF(EX(J,x).GT.EX(I,n*1))G0  TO  75 
I1=-l 
11=11*1 

SEE  IF  J SCRIBES  A 360  DECREE  ARC  ARCUNO  I 
ANG=0. 

xc=p(pi m*ii.n 

VC=P(P1 (I)*II,2) 

xp=p(P1 ( J)  » 1 I-XC 

VB=P(P1 (J),2)-VC 

Bl=xp*XP»YP*VB 

IF ( PL  .FC , 0 . )GO  TC  73 

DO  76  JlsPl CJ),P1 (J*l)-1 

XAsxp 

VAsVB 

AL*BL 

XB*P(J*l,l)-XC 

YR=P ( J* 1 » 2)“YC 

RL=XP»XP*TB*YB 

IF(PL.EC.O.)GO  TC  73 

SISE  = (X»«YB-Xe«Y*)/SORT  C*U*BL) 

AnG*AnG*ASINCSXNE> 

ANG*AB$(ANG) 

IF { ABG.IT. t )GO  TO  75 
NCONC  f“)=BCONCfM)  + l 
SUR9(NCCNC(M),P)*I 


CONTINUE 
ACCEPT  COMMANDS 
TYPE  77 

FORHATC  BRANCH*  '»*) 

ACCEPT  10.TRIC 
TPIClsTRIC 

IFfTPIC.EO.’A'JGC  TO  92 
IFCTRIC.FO.'R'IGC  TO  79 
IFfTRIC.EC.'K'IGC  TO  92 
IF(TPIC.EO.'T')TYPE  2U.NJ.NPT 
IFfTPIC.FG.'I'IGC  TO  92 
IFtTPIC.EQ.'C'IGO  TO  01 
IF  fTRIC.FO.'S'IGG  TO  110 
IFCTRIC. FO.'t'IGO  TO  139 
IF (TPIC.EQ.'T')GC  TO  76 
TYPE  70 

FORMAT ( ' AUTOMATIC, WAPP.PANUAL. INSPECT, SINGLE, CAP, EX  IT. TOTALS  * ) 
GO  TO  76 

CHANGE  NARP  ANGLE 
TYPE  00 

FORMATf  HAX.WARF  ANGLE*  ',») 

RE *rf5,*,END=76,ERPs76)ANG 
wang=C0SD ( ANG) 

GO  TO  76 

CAP 


TYPE  B2 

FOOMATC  GLOBAL  LOOP  NUMBER*'#*) 
RE  AD (5. 1 10,END  = 76,ERR*76)N 
DO  03  1*1, NLV 
J*I 

IF(LPP(I).GE.N)GC  TO  0U 
CONTINUE 


itfy 


80  IMJ.lT.?)G0  TC  86 

TYPE  BS.P(P1 flPPf J-l) ).3) 

65  FOR-ATC  7 FOR  (LEVEL-1 )»' ,E10.3) 

86  TYt>r  87,P(P| (N),3) 

67  FPR-iTf  7 FOR  (L E VEL  ) s ' , E 1 0 . 3 ) 

JF(J.Gf .KIVJGO  TC  PR 

TVPF  8P,P(P1 (LFP(J*1)),3) 

68  FORMAT!'  7 FOR  ( L E VE  L ♦ 1 ) = ' , E 1 0 . 3 ) 

66  TYPE  60 

60  F0R-4TC  ENTER  Z FOP  VERTEX  ',*) 

Rt*Pt5,*»ENDs76,fRHs76)ZV 
NJsNJ+1 

P(KJ.l)s(EX(K,?)YEX(N,l))/p. 

P(NJ.?)=fEx(N,«)4EXtN,3))/2. 

p ( n J , 3) «Z v 

pjrPJ  (NtD-PJ  (N)>1 

POsPJ (N)-l 

00  6!  I«1,IN 

JP(t,NPT)sNJ 

IP(?f RPTJsPOtl 

IP(3.NPT)sP0tI+l 

JP(0,NPT)«P 

NPT  SNPT ♦ 1 

61  CONTINUE 
GO  TO  It 

62  TYPE  63 

63  FOR-FTC  START  KITH  WHICH  LEVEL?  ’,5) 

RF40(5, U«»EH0=7t,EWRs76)LIN 
IF(LIV.EC.O)LIK=LINO 

60  IK0?=NLV-1 

IF(TPIC.FC.'I')I*'0?  = NLV 
IFfTRic.SE.'O'JGC  TO  605 
TYPE  601 

601  FORHOTC  POST-E0IT?  ',$) 

ACCEPT  10.ANS 
PEr. FALSE. 

IF(ans.EC. 'V')PE*.TRUE. 

605  00  13P  ILSLIN.INE2 

TO  = LPPt  ILyD-1 
P0*l PP ( IL )• 1 
IFCTRIC.fcE.U'JGO  TO  68 
C EMPLOY  CONNECTIVITY  ALGORITHM* 

DC  65  Jsl,8 
FLAG(J)s. FALSE. 

65  CCJ.DsO 
NL(l)sl 
NL(2)*1 
NLBsTO-eO 

KLT«LPP(IL*2)-IPP(IL6|) 

DO  67  Jsl.MLB 
DO  67  Ksl.NLT 
DO  66  Ls 1 • 3. ? 

lF(EX(PO*J,L*n.LT.FX(TO«K,LnGO  TO  67 

66  IF(FxfPO*J.L).GT.EX(TOYX,L4l))GO  TO  67 

c(j.n=cu,i)*i 

C(K*MM  )sC(«*NLPr  1 )♦  1 
C(J.C(J, 1 ) ♦ 1 )SX4NLP 
C(«66LP.C(K+NLP* 1)61 )*J 

67  CONTINUE 

IP(TRIC.EQ.,»,.ANO..NOT.PE)GO  TO  121 


c PAINT  t L*BFL  LOOPS 

98  IF(.NOT.TEKT)GC  TO  109 
00  99  1=1, u 

99  E(I)=EX(R0+1,I) 

DO  100  I=B0*2,LPP(ILt2)-l 
00  100  J=l,3,2 

1F{F(J).GT.EX(I,J))FCJ)=FX(I.J) 

100  IFfE(J*l  l.LT.EXd,  J«1  ))E(J+1  IrEXf  I,  J+I) 
OxsF  (2)-E ( 1 ) 

DYsE ( 9 ) -F (3) 

XBA9r(E(2)+E(l ))/2. 

Y5A9=(F(«)4E<3))/2. 

sf=ox 

IF(r'Y,GT.SF)SFsDY 
SFs700./SF 
I ND  = 2 

IF(TRIC.EQ.'A*)G0  TO  121 
CALL  BOX 

1F(TPIC.E0.'I')IA0=1 

101  00  107  J=1.IN0 
INOl=LPF(IL+J-l) 

IN02=LPP(IL*J)-1 
OP  107  I=IND1,!NC2 
CALL  PPVEfPl (I), 0,0) 

* 00  102  KsP1(I)Al,Pl(I*l)-l 

KKSK 

102  CALL  DPAW(KF,0,0) 

II=I-IN01*1 

CALL  «0VE(P1  CD, -57, -10) 

CALL  ALPOOF 

IF ( J.E0.2)G0  TC  105 

IFfTElC.EQ.'I'lTYPE  103,1 

103  FOR-AT ( 1 H ♦ # $ , 13) 

1FCTRIC.NE.  ' I ' ) TYPE  109,11 

109  F09-AT(1H»,S, 13, 'B' ) 

GO  TO  107 

105  TYPE  106,11 

10b  F09*'AT(lh+,3,I3,,T') 

107  CONTINUF 

CALL  mvTC(0,767) 

CALL  ALPODE 
TYPF  108, IL 

108  FOR-atC  >>  LEVEL  ',12) 

IF(TRIC.nE.'I')GO  TO  109 

RE  AO (5, 1 lu,EN0=7fc,FNRs76) ICONT 
GO  TO  138 

109  nl(1)=0 
NL(2)=0 

C SELECT  CONNECTIVITY  MANUALLY 

no  type  in 

111  FORMAT  ( 1 BOTTOM  LOOP(S)  (LOCAL):') 

GO  TO  113 

TYPE  112 

112  FOP-ATC  BOTTOM  L OOPS 5 ( GL OB AL ) ' ) 

113  REA015, 1 t«,EN0*76,EPPs76),N 

119  FOR-AT(I) 

1F(n.EL.O)GO  TO  115 
NL  ( 1 ) *NL  ( 1 ) ♦ 1 
LPSTk(1,kl(1))sN+B0 
IF(TRIC.En.'S')LPSTK(1,NL(l))eN 


•pRACT 


1CA®3* 


t 

r 


59 


115 

116 
117 


116 


119 

120 


C 

121 


122 


123 


129 


125 

126 
127 
126 


129 
C 

130 


C 

131 


GO  TO  113 

IF(NLC1).E0.0)G0  TO  119 
TYPE  116 

FOB“AT(  1 TOP  LOOP ( S ) 1 ) 

fit tO (5. 1 1U,£N0=76,ERR=76),N 

IF(N.EC.O)GO  TO  118 

Nl(2)=NL(2)+l 

LPSTK(2.ML(2))sN+T0 

IFfTPlC.f Q. 'S' )LPSTK(2,NL(2))=N 

GO  TO  117 

IF(NL(2).FQ.0)G0  TO  119 

CALL  QUAD 

TPICsTPJCl 

IF(T»IC.EQ.'A')GC  TO  138 
TYPE  120 

FORMAT!'  NEXT  LEVEL?'**) 

ACCEPT  10*  ANS 
IF(ANS.FC.'Y')GO  TO  136 
IF(.NOT.TEKT)GC  TO  109 
CALL  SOX 
GO  TO  101 

AUTOMATIC  PRANCHING  ALGORITHM 
NL(1)=1 
NL(2)=1 

IF(K'CONCCIL).FC.O)GO  TO  128 
00  127  J=1,NC0FC(IL) 
x = SURR(  J.  ID-80 
FLAG(K)s.TRUE. 

IF(NCONC(IL+1).GT.O)GO  TO  123 
!F(C(K,n.GT.nTYPE  122, IL 

FOR“AT ( ' 7ILLEAGLE  CONCENTRICITY  IN  LEVEL’,12) 

LPSTK  1 1 , 1 )=k»90 

LP6Tk(2*1)=C(K,2)aT0 

CALL  OU*D 

GO  TO  127 

DO  129  L = t ,C(K,  1 ) 

DO  129  N 2 1 * NC ONC ( IL ♦ 1 ) 
IF(SUP°(n,IL+1).E0«C(K,L*1))G0  TO  125 
CONTINUE 

LPSTK(p,1)=C(K,2)*T0- 
GO  TO  126 

LPSTk (2* 1 )=SURRCN. IL41 ) 

FLAG(SL'RR(N,  IL+1)-T0)s.TRUE. 

CALL  OUAD 
CONTINUE 

IFfNCON'CdL+l  ).EG.O)GO  70  130 
00  129  Jsl,KXONC(lL4l) 

KrSURR ( J , IL ) “80 
FL»G(K)s.TRUE. 

IF(C(K,1).GT.1)TYPE  122.il 

LPSTK ( I , 1 JSK4B0 
LpSTk (2* 1 )=CCK,2)4B0 
CALL  QUAD 

PREPARE  remaining  LOOPS  FOR  QUAD 
DO  139  J s 1 , NLB 
1FCFLAG(J))G0  TO  138 
FLAGf J)2.TPUF. 

IF(CfJ,l)-l)13e,131,132 


SIMPLE  1 ON  1 

ivsccc(J,i)4NLe,i) 


60 


IFdV.GT.DGO  TO  132 
FLAGfCU, 2)1®. TRUE. 

K‘L  ( 1 ) s 1 
MC?)  = 1 

LPSTM?.  1)=C(J,2)4T0-NL8 
LPSTk(1,1)=C(C(J,2),2)4B0 
CALL  QUAD 

go  to  lie 

C DETERMINE  EXTENT  OF  BRANCHING 

152  SP=1 
NBR=1 

BRANCH  ( 1 ) s J 
STACK(1)=J 

153  DO  13a  I=1,CCSTACK(SP),1> 
N0DE=CfSTAf:*(SP),l4l) 

if(flag(nodE))go  to  13a 
FLAG(NOCE)  = .TRL'E. 

NRR=NRR4l 
BRANCH(NPR)=NODE 
IF(CfNOCE.l).EC.l)CO  TC  13a 
SP=SP+1 

STACK(SP)=NODE 
GO  TO  133 
13a  CONTINUE 

IFCSP.EC.t )G0  TO  135 

SP=SP-1 

GO  TO  133 

135  NL(l)sO 
NL(2>=0 

00  137  1=1, NRR 

IF(PRANCH(I).LE.NLB)GO  TO  136 
NL(2)=Ht2)4l 

LPSTt«(?,NL(2n  = BfiANCN(I)4T0-NLB 
GO  TC  137 

136  *'L  ( 1 ) *NL  ( 1 ) ♦ 1 

LpSTK(1,NL(l))  = BRANCH(I)4B0  <• 

137  CONTINUE 
CALL  QUAD 

136  CONTINUE 

IF(tric.EO.'I')GO  TO  76 
LIND=NLV*1 
IF (OAT AF )G0  TO  6 
C OUTPUT 
139  TYPE  1 ao 

iao  FOR-ATC  OONE?  '.*) 

ACCEPT  10.ANS 

IF(ANS.NE. »V* )G0  TO  6 

NPTbNPT-1 

NPsNPMAX-1 

TYPE  2a,Nj,NPT 

IPLC8.IPLI.NIPL(IPLI))«NPT 

TYPE  1 a 1 

la  1 FORMAT  C ' OUTPUT  FILENAME?  ',S) 

ACCEPT  a. CHAPE 
OPEN(UMT  = ??,FILE«ONAPE) 

DO  la?  1=1, NP 
NP I C ( I ) = 0 

OC  la?  J*1,NIPL (I) 

la?  NPiC(n=NPICCI)*IPL(2,I,J)-IPLf  1,1*  J)*l 
NPlTE(22,la6)NP,NJ,NPT 


61 


Nr  0 

DO  103  IU.NP 
J=NM 

10 3 npla(2.I)sN 

NRITt (22. 108) r (NFLA(I,J),Iri,2),Jsl,NP) 

NBJTE (22,107 )((P{ I. JJ.Jr  1,5), Irl.NJ) 

nRITE(22#106)((<(IPCI,J)#Is1#0),Js:pl(1,N,H),IPL(2,N,M>) 

1 ,rsl,MPl(N)),N«l,NP) 

TYPE  100 

too  for-atc  special  function  file?'#*) 

ACCEPT  10.ANS 

if(ans.ke.'t')go  to  ioe 

TYPF  105 

105  FOR-ATC'  FILENONE  OF  S.F.  FILE:  ',*) 

ACCFPT  O.Sna-E 

OPf N fUNITe?J,P ILEsSNONE) 

n»ITE(?3. 107) (PCI, 3).I*1#NJ) 

10b  FORMAT ( 20 1 0) 

107  F OR- AT (fc£ 12 .5) 

108  STOP 
END 

subroutine  quad 

LOGICAL  Ft  AG1 ,FL*G2. FLAG (5). CLOSE. IPF (1 0 ) , ORAhF .MAP.eo , TEKT.PE # PEN 
C0““0N/0/M  (2)#LFSTK(2.5)»NPT»N,J#02»B0»TEKT»PE. 

1 E*(100,o),Pl (100).TBIC.IpLf2,5,10),MPL(5),IPLI#NP*AX,IP(0,2000) 

CO‘<*CN/NOOFy’B2(2, 1 100,3) 

CfM-CN/TM/XBAfi,  ybab.sf  ,p(2000#3)#  IYD2 

INTFGFB  TEHPC200).NV(2),NEKC2,5), V (2, 200) , F IRST, PI ,PDUB C 12) 

REAL  TOAKSC2,U),1mAP12,2),S12) 

INTEGER  CRDfcR(5,5).CP(S.S),OC5),FIRSTl ,F IHST2 .CL.CT 
PEAL  0(5,5).CTf5).E(0) 

INTEGER  NCP(2),CVC2, 150),N1 (2),N2C2) 

DATA  IVIt>,*AP,IND6/l,  .TRUE. ,2/ 

NPTlsNPT 

NPT2rNPT 

NNEWSO 

NTEMPsO 

XBAPTrXBAR 

YRARTs YB AR 

NTHAXst 

NB-A  X e 1 

NLINE=767 

PENS. FALSE. 

SFTsSF 


DRAUFs. FALSE. 

DO  OB  Nr  1 # 2 
OCl)*t 

IF(NL(N).E0.1 )G0  TO  2R 
IF(TRIC.NE.'a*.OR.PE)GC  TO  12 
C AUTOMATICS 

C find  CLOSEST  OISTANCE  BETWEEN  ALL  LOOPS 
DO  7 Irl,NL(*)-l 
FlPSTleP)  (LPST«(F,n) 
l»ST 1*P1 (LpSTK{M, I ) ♦ 1 ) • 1 
DO  7 JsI«l,NL(P) 

FIPST2rRl (LPSTF(N.J)) 

L»ST2spl (LPSTK(w, J)*l)-1 
CPU.  J)sFIRSTl 
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t 


2 


4 


5 

6 
7 
C 


8 


4 


10 

11 

c 

12 


CP(j.nsFinsT2 
0TST=1E35 
Fl*Gl=. FALSE. 

Fl*G?s.1BOF. 

00  3 K = F I 95T 1 » L AST  1 

TDISlsCPt*. 1 )-P(CP(J,  1),  1))**24(P(K,2)-P(CP(J»I).2))**2 
IFdOIST.GE.DISTlGO  TO  3 
DO  2 L=1,NL(h) 

1F(k.E0.CB(I,L))G0  TO  3 

PIST=T01ST 

FL*G2=. FALSE. 

CP(I.J)=K 

CONTINUE 

IF<FHG2.AN0.FLAG1)G0  TO  6 
FLAG1=.TRUE. 

DO  S K=FI«ST2,LAST2 

T01ST  = tP{K,l).P(CP(I,J),m**2*(P(K,2)-P(CP(I,J),2))=»2 
IF(TOISI.GE.CIST)GO  TO  5 
DO  u l=l,Nl(P) 

1F(k.E0.CP(J,L))G0  TO  5 

DISTrTDIST 

FLAG1  = . false. 

CP(J,n=N 

CONTINUE 

IF(FL*G1  .AN0.FL*G2)G0  TO  b 

GO  TO  1 

DCJ,I)=DIST 

D(I,J)=DIST 

PICA  The  PROPER  loop  sequence 
DO  9 I=1,NL(H> 

J1*I 

DT  1 1 ) = 0 . 

00  8 J=1 ,5 
FLAGUJs.  false. 

DISTsO. 

00  9 J=1,NL(H) 

FLAGf JlJs.TBUE. 

DT ( I ) =0T ( Il+DIST 
ORDER ( I , J ) = J 1 
J2  = J 1 
DIST=1E35 
DO  9 K=1,NL(H) 

I F (F  L *G ( K ) ) GO  TO  9 
IF(OIST.LE.OCJ2,K))GO  TO  9 
JlsK 

DIST=D(J2,«) 

CONTINUE 

PISTs1E35 

DO  10  I=!»NLtH) 

IFfOIST.LE.DT(I))GO  TO  10 
DISTeDTCI) 

IT  = I 

CONTINUE 

DO  11  Jsl.KLOO 

0(J)=0RCER(IT,J) 

OISLPAF  LOOPS 

IF(TRIC.EG.,A,.ANO..NCT.PE)GO  TO  26 
1F(.N0T.T£KT)GC  TO  17 
CALL  BO* 

DO  16  Ks 1 , 2 


W'8 
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t 


00  16  Jsl.MOO 

1 NO  1 sP 1 (LPSTK(K,J)) 

IM0?rPl CLPSTKfK, J)+l)-l 

CALI  **OVE(IM01,0,0) 

CO  13  L*IKDltl,IKD2 

ll  = L 

13  CALL  DRAl»(LL»0»0) 

IFC*.NE.K)G0  TO  16 
CALL  ALPODE 
TYPE  1 A>  r J 

10  FORMAT (IP*, 3,  *L0CP*»I1) 

00  16  L=IND1+2,IND2,2 
Ll=L 

CALL  CPCSS (LL , 3 ) 

CALL  POVECLL»-57,-10) 

CALL  ALPODE 
TYPE  15,1 

15  FOPPAT  Cl  *♦ , J. 10) 

16  CONTINUF 

CALL  MVT0(0,767) 

CALL  ALPOOE 

17  TYPE  18 

18  FOPPATt*  ENTER  LCOP  SECUENCE') 

N = 0 

DO  10  L*1,NL(M) 

NsN«l 

1R  READ(5.«.ENOsl37.ERR=137)0(L) 

20  NL(P)»N 
TFCN.EO.DGO  TC  26 
00  23  L=1,NLCM)-1 
TYPE  21.0CL) 

21  FORMAT  ( ’ L'.Il,'  P?',») 

ACCEPT  •,CP(0(l),0(L+l)) 

TYPF  22,0(LO> 

22  F0O“AT(lM  + ,$,  • TC  L ' , 1 1 , 1 P?'.S) 

23  ACCEPT  *,CP(0(L*1),0CL)) 

TYPE  20 

2o  for-atc  changes? ' , $) 

ACCEPT  25, ANS 

25  FOR-AT(Al) 

1F(AKiS.EC.«Y')G0  TO  17 

C INTERPOLATE  NE*  POINTS 

26  IFCNICm).EC.1)G0  TO  2R 
IsLPSTK(2,l) 

261  JsLPSTK (1,1) 

262  llsPl(I) 

263  JlsPl(J) 

2 60  71sP(Il,3) 

265  ZJs°(Jl,3) 

266  Z = (7WJ)/2. 

00  26  I«1,NL<H)-1 

NJ*NJ*1 

NEP(P,I)SNJ 

DO  27  J:  1 , 2 

KsJ 

270  I?*CP(0(I)»0CI*1)) 

271  !3*CP(0(!*1),0(I)) 

272  TJ?aP(I2, J) 

273  TJ3*P(I3,J) 

270  TJ«CTJ2*TJ3)/2. 


jfr' 
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27  P(KJ,J)sTj 

28  P(NJ,3)sZ 
C COMBINE  LOOPS 

29  F l AG  1 = . TRUE . 

FLAG2=. FALSE. 

N V ( M ) bO 

I  NO l * I 
INC2*NL(M) 

INDJbI 

JO  00  03  1 1 1 N01 , I NO? « IND3 

JI  FIRST=P1(LPSTK(M,0(I)))-1 

IF(M(»*).E0.1)G0  TO  JR 
CL=CP(0(I),0CI«n)-FIRST 
IF(I.FC.l.AND.INC3.E0.-nG0  TO  38 
IF ( I .EQ. 1 )G0  TC  32 
CTbCP(0CI)»0(I-1J)-FIRST 
IF (I.FO.Ino?.anD.IND3.EQ.I)GO  TO  35 
IF(.noT.FLAG1)G0  TO  33 

32  FLACls. FALSE. 

FLAG2=.TRUE. 

N = CL 

IFCIMJJJNaCT 
GO  TO  OO 

33  IF(IN03)36,3O 

30  IF(CL.GT.CT)GO  TO  37 

F L AG  1 = .TRUE . 

35  F IRSTrFIRST+CT-l 

NsPJ (LPSTK(M,0(I))*l)-FIRST-2 
GO  TO  00 

36  IF(CL.GT.CT)GO  TC  38 

F IRSTrF IRST+CL-1 

N=C  T— CL ♦ 1 
FL»G2=.TRUE. 

GO  TO  OO 

37  FIRST=FIRST+CT-1 
FL*G2».TRUE. 

NsCL-CT+i 

GO  TO  00 

38  FIPST=FIPST*CL-1 
FL  *G 1 = .TRUE • 

39  I!H<»-  = LPSTK(M,0Umi 
N=P| fI0LM)-FIRST-2 

C ALL  That  hOSK  FOR  THIS  LITTLE  LOOPI 

00  IF  C I . EC. IN02JFLAG1S. FALSE. 
IF(N.EO.O)nb1 

00  0|  J=1,N 
NV(M)rNV(M)M 

01  V (MfNV(Ml )pF I RS  T* J 
IF(.NOT.FLAG2)GO  TO  02 
N=NEN(M,I) 

IFnNr'3.E0.«l)NsNEN(M#I*l) 

NVCM)=NV(M)»1 
V(m,nv(m))bN 
FLAG2=. FALSE. 

02  IF (FLAG! )G0  TO  31 

03  CONTINUE 
1F(NL(m).E0.1)G0  TO  00 
IF(IND3.E3.-1)G0  TO  Ofl 
FLAG  1 5 .TRUE • 

1*033-1 
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65 


i 


as 

C 


45 

46 

47 

48 

c 

49 


50 


51 

52 
55 

54 


55 

56 

57 


INOIsnl(h) 

IND2*1 
GO  TO  30 
CONTINUE 

SEGMENT  ENDPOINTS 
NV(M)sNV(M)*| 

GO  TO  46 
DO  45  1=1.2 

IF(P(V(M,n,I).NE.P(V(P,NV(M)),I))GO  TO  47 
V(M,NVp>))sV(N,  1) 

GO  TO  4e 

V(M,NVfM))sV(M,NV(M)-l)*l 

CONTINUE 

IE(TPIC.EQ.U')TRIs**» 
IF(TRlC.EC.,A'.AND..NOT.PE)GO  TO  111 
DETERMINE  EXTREMES 
DO  49  1=1,4 
EC1)=EX(IPSTK(1,1),I) 

DO  50  M=l,2 
00  50  J=1,NL(M) 

00  50  1=1, 3, 2 
ET1=EX(LPSTK(M,J),I) 

ET?=E*(LPSTK(M,J),Itl) 

IF(ETl.lT.E(I))E(n=ETl 

IP(ET2.r.T.E(I«l))E(l4|)sET2 

0X«E(2)-E(1) 

0Y=E(4)-EC3) 

OYT  sOVtDZ 
SF?=DX 

IF(dyt,GT.$F2)SF2=DYT 

SF?s700./SP2 

SFlsDX 

lF(DY.GT.SFnSFl=DY 
SF 1S700./SF1 
XRAPs(E(1)+E(2))/2, 

YRARsCE(3)*E(a))/2. 

IFnviEN.ED.DGO  TO  52 
SFsSF? 

IYD2=390.*OZ/0YT 
GO  TO  53 
SFsSFl 

1 Y02  = 0 

IFf.NOT.TEKDGO  TO  63 
CALL  BOX 
00  54  M= 1 , ? 

C*tL  CRCSSfV (m, 1 ) ,5) 

00  54  J=2#NV(M) 

C»LL  DRAN(V(M, J),0,0) 

IFCIVIEx.EO.nGO  TO  57 

1 Y*  1 0 

00  56  Msl,2 

IXs$F«(P(v(M, D, 11-XBAR )*634 
CALL  MVTO(IX.IY) 

DO  55  J=2,NVfM) 
IXrSF«(P(v(^,J),l)-XBAR)A634 
CALL  VCTO(IX.IY) 

IYrSF«0X*10 


IFCNPTI.FO.NPTJGC  TO  60 

IN0=NPTJ*1 

DO  56  IsIND.NPT 


c*ll  movecipci,i),o,o) 

58  CALL  0RANltPC2,I),0,0) 
IFCIVIEu.FC.l)GO  TO  60 
IYTs<i»IyD2/39 

I YBs35* I YD2/39 
DO  59  IsInC.NPT 
ix*SF*(P(iPci,i>.n-xeAf»)*65o 
call  mvTOCIX.IYP) 

IxsSF.CPHPC2,I),  1)-XPAR)*630 

59  CALL  VCTOCIX.IYT) 

60  IF CTRJC.EQ. ' A ' .ANO..NOT.PFN)GO  TO  111 

CALL  mOVFCVU, 11,-10,-10) 

CALL  ALMOOE 
TYPE  600 

600  FORMATOH*,*, 'e') 

CALL  KOVECVC2.1),-10,-10) 

TYPE  601 

601  FORmaTCIH*,*, 'T') 

DO  62  Ms  1 f 2 
J=NV(M)-2 

00  62  IsM,J,INC6 

CALL  MOVE(V(M,I). -«5,-10) 

CALL  ALMODE 
TYPE  61,1 

61  F0RmaTC1M*,J,I3) 

62  CONTINUE 

CALL  MVTO ( 0,767 ) 

NLINE=767 
CALI  ALMCDE 

63  1P(ORAKF)GO  TO  90 
C ACCEPT  COMMANDS 


60 

65 


66 


C 

67 
C 

68 
C 

69 


TYPE  65 

NL I N£  sNL I NE-20 


FQRmatC  TRIANGULATE>',S) 


ACCFPT  136 

, TR I 

IFCTRJ.EC. 

' A ' )G0 

TO 

111 

IFCTRI.EC. 

' V'  )G0 

TO 

72 

1FCTRI.EC. 

•R')GO 

TO 

90 

1FCTPI.FC. 

’C')GO 

TO 

137 

IFCTOI.EC. 

' T ' )G0 

TO 

73 

1FCTPI.FC. 

• I ' ) GO 

TO 

9il 

1FCTPI.EC. 

'O' )G0 

TO 

69 

IFCTRJ.EC. 

•P'JGO 

TO 

75 

IFCTPI.EC. 

•E'JGO 

TO 

71 

IFCTRJ.EC. 

•O' )G0 

TO 

83 

IF  CTRI.FC. 

•M' )G0 

TO 

67 

IFCTRI.EC. 

' U ' ) GO 

TO 

68 

TYPE  66 


NLINEsNLINE-80 

FOPMATC  auto, INTERACTIVE, ',/,  ' RE NUMbE R , DE NS I T Y , • , /, • CONT, 
1 TOTALS',/,'  PART, ERASE, VIEN, MAP, UNMAP, ONE') 

GO  TO  60 
SET  MAP  FLAG 

MAPs. TRUE. 

GO  TO  69 

RELEASE  map  FLAG 
“APs. FALSE. 

GO  TO  60 

SET  NU-BER  DENSITY 
TYPE  70 


V\SP 


67 


t 


70 

C 

7 1 


C 

72 


C 

73 

74 

75 
C 

76 


77 

78 

79 

80 


81 

82 


C 

83 


84 

85 


8b 

87 

88 


FORMAT!'  1 ND3* ' # * ) 

ACCEPT  *» IND6 
CO  TO  53 

ERASE 

NPT«NPT1 

NTMAXsl 

NBMAXa] 

ORAwFs. FALSE. 

GO  TO  53 

CHANCE  NU'-BER  OF  ORTHOGRAPHIC  VIEWS 


1*1  VIEW 

IF(I.E0.2)IVIEW»1 

IF(I.EC.l)IVIEWs2 

GO  TO  51 
TOTALS 

TVPF  74,NJ,NPT 

FORMAT!'  ',14,'  NODES • ',14,'  ELEMENTS') 

GO  TO  60 
TYPE  76,IPLI 
CHANGE  PART  NUMBER 

FORMAT!  ' CHANGE  FART  FROM  Ml,'  TO  ',*) 

ACCEPT  *,N 

IF(N.LE.NPMAX.ANC.NPT.GT.lPL(l»IPLlrNlPL(IPLI)))GO  TO  80 
IF(N.GT.NPMAX)G0  TO  78 
TYPE  77 

FORMAT!'  CANT  CHANGE  PART  NUMBER  YET') 

GO  TO  64 
TYPE  79  # NPMA  X 

FORMAT!'  MUST  PE  LESS  THAN  ',12) 

GO  TO  75 

IFIIPLI.FO.N)GO  TO  6« 


IsNPT-l 

IF!IPL(1.IPLI,MPL!IPLI)),LT.I)G0  TO  81 

NIpl!IPII)=nIPL!IPLI)-1 

GO  TO  82 

IPL!2,IPLI,MPL!IPLI))*I 
NJPL  !N)sMpl  !n)41 
IPL!l,N,NIPL!N'))sNPT 
IF !N.EQ.NPmaX)NPHAX=NPMAX+1 
IPLI=N 
GO  TO  64 
ONE  AT  A TIME 

00  88  1*1,4 
NL INEsNL INF-20 
TYPF  84,1 

FORMAT!'  POINT  ',11,'.  T/0  » NO. ' » S ) 
ACCEPT  85, T, N 
F ORMAT  ! A 1 , 1 ) 

IF!T.EO.'T')GO  TO  86 
I F ( T . E 0 . ' B ' ) GO  TO  87 
IF!I.NE.4)G0  TC  64 
IP!4,N»T)*0 
GO  TO  88 
IP(I,NPT)sV!2,N) 

GO  TO  88 

IP!I,NPT)sV!l,N) 

CONTINUF 

CALL  movE!IP!1,NPT),0,0) 

CALL  OR*w!IP!P,npt),0,0) 

CALL  DPAW!IP!3,NPT),0,0) 


IFUPCO.npt).EC.0)CAU  DRAWCIPCO,NPT),0,0) 

*9  NPTsNPTa  1 

CALL  HVTOCO, KLINE) 

CALL  ALPCOE 
GO  TO  60 
C RENUHRFR 

90  TYPE  9 l 

91  FORRATC  NE h *1  NODE',/*'  HOP  CONTOUR  )',*) 
NLJNE=NLINE-«0 

ORA*F  = . FALSE. 

ACCEPT  *»N5 
IF (N5.LT . 1 )G0  TO  90 
GO  TO  113 

90  NTMTsNTPAX 

NBMTsNRPAX 
TYPE  95,NTRAX,NV(2) 

NL1 NE  sNL INE-80 

95  F OR*  A T C ' TOP: ' , 13. 13) 
READC5,96,EN0=6«,ERR=60),M C2),N2C2) 

96  F0R“ATC2I) 

IF(N1 C2).E0.0)G0  TO  137 
TYPE  97,NB“AX,NVC1) 

97  FOR“A T ( ' POTTO*:', 13. '-*.13) 
RF»0C5,96,EN0s60,ERR=fc0) ,N1 ( 1 ),N2C1) 

NPT2=NPT 

98  IFtTRI.NE. 'A* )G0  TO  99 
Nl(l)=l 

N 1 C2)  = l 
N2  C 1 ) = NV  1 1 ) 

N2C2)=NVC2) 

99  NT0EL»1 

IF  CM (2) .GT.N2(2))NTOELs-l 
NB0EL=1 

IF  { N 1 ( 1 ) .GT.N2C 1 ))NRDEL=-1 
IF  (.►. OT.*AP)GO  TC  108 
00  100  1=1,2 

100  !F(M  U).EO.N2(I))G0  TO  110 
00  102  1=1,2 

I N D 1 =N 1 (1) 

IND2=N?(I) 

IF  (U01  .LT.IN02)G0  TO  101 
I ND 1 = 1 N02 
IND2SN1 (I) 

101  TRAN SCI, 1)=P(V(I,IN01),1) 
TRANSCI,2)=TPAKS(I,1) 
TRANS(I,3)=P(VCI,IND1),2) 
TRANSCI,R)=TRANS(I,3) 

00  102  J=IN01A1,IN02 
PXsPC  VCI , J), 1 ) 

1FCPX.LT. TRANS(I,1))TRARS(I,1  ) = PX 
IF CPX.GT. TRANS C 1,2) ) TRANS  Cl, 2)*PX 
PY  = P C V C I , J ) , 2 ) 

IF(PV.LT. TRAN SCI, 3))TPANSCI,3)=PY 

102  IF  CPV.GT.TRANSC  1,0))  TRAN'S  CI,«)=PY 

103  on  too  1=1,2 

T«APC 1 , I ) = CTRAN  SC  I ,?) *TR AWS  C I , 
too  T*AP(2, 1)  = CTPANSCI.3)*TPAN'SC  I, 

DO  105  1=1,2 

IFCTRANSCI,2).E0.TPANSCI, 1 ))R0 
105  IFCTRANSCI,3).E0.TRANSCI,«))G0 


I / / * C • 
0))/2. 

to  toe 

TO  10B 


& 


C MAP 


00  107  1*1,2 

Sms>nO./fTfiA*.5(I.2)-TRANS(I»n) 
S(2)=100./(TRANS(I,R)-TRANS(I, J)) 


INSKPOEL 

IFd.E0.2>IN=NTDEL 
DO  107  J*Md),N2(I),IN 
00  106  K*l,2 

106  P?tl, J.K)aS(K)*(P{V(I,J).K)-TKAP(K#I)) 

107  P2(I. J,3)*P(V(I,J),3) 

GO  TO  110 

C DON'T  “AP 

108  DO  109  1*1,2 
INsNPDEL 

IPd.EO. 2JINSNTDEU 
00  109  J = N1  (I),N2tn,IN 
DC  109  K*l,3 

10O  P2d,J,K)sP(V(I,J),K) 

110  NT  SN 1(2) 

N6  = N 1(1) 

NTWAX=N2(2) 

NP*>AX  = N2(1) 

GO  TO  117 

111  D1S1E35 
*l=P(Vd,l),l> 

TlsP(V(l,l),2) 

DO  112  I=1,NVC2) 

02=(PfV(2,I),  l)-Xl)**2*(P(V(2,n,2)-Vl)«*2 
1F(C2.GF.D1)G0  TO  112 
01*02 
N5=I 

112  CONTINUE 
C RE-0RDE9 

113  IF(N5.rr.l)G0  TO  116 
DO  1 1 R J*1,NVC2) 

110  T f MP (J)*V(2»J) 

DO  115  J*1,NV(2) 


JR*J*N5-1 

IF (JR.GT . N V ( 2 ) 1 JR* JP-NV ( 2 ) ♦ 1 

115  V(2,J)sTE**P(JR) 

116  IFITRI.FG.'R'IGO  TO  53 
GO  TO  Re 

C triangulate  between  limits 

117  IFlAGsO 

118  7°  f 1 » NPT )*V ( 2, NT ) 

IP(2.NPT)SV(1,KB) 

IP(0,NPT)sO 

PX*P2(?,NT,1 )-P2(l,NB,l) 

pysP2(2,m,?)-F2(1,nB,2) 
R7*p?(?,nT,3)-F2(1,NB,3) 
IF(NT.EC.M“AX)GC  TO  125 
IF(NB.EG.N8MAX)GC  TO  120 
01=0. 

D2*0. 

DO  119  1*1,3 

01*01 ♦<P2(2,nt, I )-P?( 1 ,NP+NRDEL, I) )**2 

119  i.2*D2»CP2(2,NT4\TDEL,I)-P2d.NB,I))**2 
C SELECT  SHORTEST  DIAGONAL 

IF ( D I .LT.021G0  TC  125 
C TOP  NODE  IS  CLOSEST 


4?/ 


& 


A 


r 

f 


f 


120 

121 


122 

123 

12« 


C 

125 

126 


127 

128 

129 


130 
C 

131 

132 


kTsNT*NTDEL 

JP(3.KPT)aV(2,KT) 

KPT  sNPT*  1 

S*  = P?(2»M,n-P2(l,NP,l) 
SYsP?(?,KT,2)-P2(1,nB,2) 
S7sP2C2,M,3)-P2tl.NB,3) 

XlrBY»S2-«2*SV 
V 1 =P7«SX-RX«SZ 
21rBX*SY-«Y*SX 
»RF*lsfXI»Xl*Yl*Vl*21*21) 

1F(AFEA1.NE.0)G0  TO  122 

NPT=KPT-1 

GO  TO  118 

1F(1FL*G.E0.2)G0  TO  12« 

JFL*G=1 
GO  TO  118 

A\GT=(Xl*X2*Yl*Y2tZl*Z2)/SQBT(AREAl*ABEA2) 

IFf AKGT.LT.t.AKG)GO  TO  123 

IFIAGsO 

10fu,NPT-2)sV(2,KT) 

K=T=NPT»1 
GO  TO  118 

BOTTO-  MOPE  IS  CLOSEST 

IF  (►'B.EC.KB“AX)GC  TO  130 
NBsk  B + NBDEL 
IP(3*KpT)sV(l#KB) 

NPTrKPT+1 

SxsP2(lfKB,l)-F2(l,KB-KB0EL.l) 

SY=P2CI ,Ke.2)-F2(l,KB-KBD£Lr2) 

S2=P2(l,Ke,3)-P2(l.wB-KBDEL*3) 

X2rPY*S2-02*SY 

y?=p?*SX-bx«SZ 

72=P>»SY-RY«SX 

ARE A2=(X2*X2+Y2»Y2»22«22) 

IF(apf A2.NE.0.JGG  TO  127 

npTskPT-1 

GO  TO  118 

IF  (IFLAG.EO.DGO  TO  129 

IFLAG=2 

GO  TO  118 

akgts  (XI  »X2*Y1«Y2*7  W2)/SCRT(AREA1*A»EA2) 
IF(AKGT.LT>ANG)CO  TO  128 
1P(<j,*.PT-2)  = IP(  l.NPT-1) 
Ip(3**»PT-2)aIP(3,NPT-l  ) 

NPTrKPT-1 
JFLAGsO 
GO  TO  118 

IF(T.  iC.EC.,A'.AKO..NOT.PE)GO  TO  137 
P A I KT  ME  to  DIAGONALS 

1F(.kOT.TFKT)G0  TO  13<l 
DO  132  laNPT2,K'PT 
CAIL  P0VE(IP(1,I),0,0) 

C*LL  99AW(IP(2,I>,0.0) 


IF ( I V IFw.EG. 1 )GO  TO  1 39 

I VT  sil  • I TT2/39 

lYf=35»IYT2/39 

00  li3  laNPT2,KPT 

lxsSF*(P(IP(l,I),l)-XMR)t63« 

CALL  PVTCCIX.IYB) 

Ix=sf»(P(1P(2.1).1)“XBAR)»63R 


135  CALL  VCT0(1X,1YT) 

13a  NLINEsnlinE-100 

CALL  HVTO ( 0*  NLINE ) 

C*LL  ALMODE 
TYPE  135 

135  FOR“ATf  CHANGES?*,*) 

• ACCEPT  136, ANS 

DRAhFs.TRUE. 

IF(n'LIne.LT.200)GO  TO  53 

136  FOPHATCA1) 

IF(AHS.EO.*Y*)CO  to  138 
IF(TPI.EG.»A')GO  TO  137 

1F(NTMAX,EQ.NV(2).AND.NBMAX.EO.NVC1))GO  TO  137 
GO  TO  aa 

138  nthaxsntht 

NBHAXsNBMT 

NPTsK'PT? 

IF(TPI.ER.*A*)CRA-F=.FALSE. 

IF(T»IC.AE.'A*>GO  TO  53 
PEN*. TRUE. 

TrPE  139 

139  FORMAT ( 1 CHANGF  CONNECTIVITY?  *,*) 

ACCEPT  136, ANS 
IF(*ns.NE.’Y*)GO  TO  53 

TRICs'H* 

C RESTORE  PARAMETERS  BEFORE  RETURNING 

137  YBARsYBART 
XBARsxBART 
1 YD2»0 
SFsSFT 
RETURN 
END 


SUBROUTINE  BOX 

CALL  CLMOA 

CALL  HVTO (290,779) 

CALL  VC YC (1 023 « 779) 

CALL  VCT0(1023,0) 

CALL  VCTC(2aa,0) 

C*LL  vCT0(2ao,779) 

PETURN 

END 

SUBROUTINE  CROS  S ( 1 , 1 SI ZE  ) 
COhhCn/TFK/XPaR , YBAH,SF,P(2000, 3), IYD2 
CALL  “OVE (I,-ISIZE,0) 

CALL  ORAh(I.ISIZE.O) 

CALL  “OVF  fl.O.ISIZt) 

CALL  DRan(I,0,-ISIZE) 
call  MOVF  (1,0,0) 

RETURN 

END 

SUBROUTINE  HOVE (1,1X0, IYD) 
COHHON/TFK/XBAB, YR*«,SF,P(2000,3), IYD2 
i x=SF«(pn,n-xeAB)  *630*1x0 
IY=SF*(F(1,2)-YBAR)*390*IYO*IYD2 
call  HVTOCIX.IY) 

RETURN 

END 

SUBROUTINE  ORAN (I, 1X0, IYO) 
COHH0N/1EK/XBAR, YBAR,SF ,P(2000, 3) , IY02 


72 


IX  = SF*fPCI,l)-)iBAR)*65«UlXD 
IVsSF  •(?(!, ?)-YB*R)  + 3<»0'UYO»I  YD? 
CALL  VCTO(IXflY) 

RETURN 

END 


APPENDIX  B 
USER  DOCUMENTATION 


Throughtout  this  documentation,  cues  typed  by  the  computer 
are  mderlined.  The  program  commences  by  asking: 

BRAIN  DATA? 

User  replies  "yes",  or  "ro"  (brain  data  is  handled  slightly 
differently  than  user  generated  data).  Next,  the  program  requests 
FILE  NAME  OF  INPUT  DATA? 

Whereupon  the  user  types  the  file  name.  The  brain  data  files 
are  called  Bn. DAT  where  n is  an  integer  from  1 to  22.  Now  the 
program  requests  a command  by  giving  the  prompt: 

READ> 


This  prompt  ;?s  so  named  because  it  is  here  that  parameters 
are  set  preparatory  to  reading  the  data  files.  READ>  responds  to  the 
following  commands,  of  which  only  the  first  letter  is  required. 


PARAMETERS  - 

Chanqe  node  elimination  parameters  S . , S 
.a  min’  max’ 

and  0 . . 
min 

TOTALS  - 

Causes  the  current  number  of  nodes  and  panels 
to  be  output  on  the  terminal. 

SCALE  - 

Enables  change  of  scale  factor.  Default  is 
.0001  for  brain  data  and  1 for  non-brain  data. 

CLOSE  - 

Computer  responds,  CLOSE  ALL  LOOPS?  If  yes, 
every  contour  line  will  be  treated  as  a closed 
loop.  Otherwise,  only  contour  lines  whose 
first  and  last  nodes  agree  will  be  considered 
closed.  Default  is  "yes". 

DEVICE  - 

Computer  responds:  TEKTRONIX  SCOPE?  "Yes" 
enables  graphical  output.  "No"  prohibits  it. 
Default  is  "yes". 

CLOCKWISE  - 

Computer  responds:  CLOCKWISE  ORDERING?  "Yes" 
forces  all  loops  to  run  clockwise.  "No"  forces 
counterclockwise  ordering.  Default  is  "yes". 
Clockwise  loops  result  in  counterclockwise 
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EXIT  - 
LEVEL  - 


BRANCH  - 

MANUAL  - 


definition  of  visible  panels.  This  understand- 
ing is  useful  in  using  the  Poor  Man's  hidden 
surface  elimination  in  MOVIE. BYU. 

Exits  from  the  program. 

Since  there  are  no  convenient  default  values 
for  the  parameters  set  by  the  LEVEL  command, 
this  is  the  command  which  actually  initiates 
the  reading  of  the  disk  file  - after  the 
following  two  questions  are  answered: 

Z SPACING=  . User  responds  with  an  integer  "n" 
which  indicates  that  every  "nth"  contour  level 
is  to  be  read  in  from  disk.  The  computer 
then  asks:  LEVEL  RANGE. -.The  user  enters  two 
integers,  j and  k.  Thereupon,  the  algorithm 
procedes  to  read  every  nth  contour  level 
from  j to  k inclusive,  imposing  the  node 
elimination  parameters.  When  array  P has  been 
loaded,  ready  for  triangulation,  the  program 
issues  the  BRANCH>  prompt. 

Proceeds  to  the  BRANCH>  prompt  without  reading 
from  d i sk. 

Same  as  LEVEL,  except  loops  from  the  brain  data 
can  be  constructed  interactively  from  their 
constituent  segments.  This  command  is  never 
needed  except  in  cases  where  closed  and  non- 
closed  contour  lines  coincide  on  the  same  level. 


The  commands  of  PARAMETERS,  TOTALS,  SCALE,  KLOCKWISE,  CLOSE, 
and  DEVICE  return  to  the  READ>  prompt.  LEVEL,  BRANCH  and  MANUAL 
proceed  to  the  next  prompt: 

BRANCH> 


BRANCH^  responds  to  the  following  commands. 


TOTALS  - 


Same  as  above. 


EXIT  - 


Same  as  above. 


AUTOMATIC  - Computer  responds 

START  WITH  WHICH  LEVEL? 

User  responds  with  an  integer,  n. 

The  computer  asks, 

POST-EDIT? 

User  responds  "yes"  or  "no",  the  consequence  of 
which  will  be  explained  momentarily.  The 
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WARP  - 


INSPECT  - 


MANUAL  - 


branching  algorithm  is  now  invoked  which 
determines  connectivity  on  the  basis  of 
window  overlap.  Branching  and  triangulation 
are  also  handled  automatically.  If  POST- 
EDlTing  was  refused,  the  entire  data  file 
is  thus  triangulated  - with  no  gaurantees. 

If  POST-EDITing  was  requested,  each  triangu- 
lated contour  pair  is  displayed  for  approval 
as  the  computer  asks, 

CHANGES? 

Two  things  may  be  in  error  - the  triangulation 
itself,  or  the  choice  of  loops  involved  in 
branching.  If  either  of  these  errors 
necessitate  change,  answer  "yes".  The 
computer  then  asks, 

CHANGE  CONNECTIVITY? 

If  the  branching  was  at  fault,  answer  "yes". 
This  has  the  effect  of  changing  from 
AUTOMATIC  to  MANUAL  mode,  and  the  user  is 
free  to  manually  control  the  branching.  If 
only  the  triangulation  was  bad,  answer  "no". 
This  moves  the  program  counter  to  the 
TR I ANGULATE>  prompt,  so  the  user  can  guide 
triangulation  manually. 

If  no  changes  at  all  are  needed,  a carriage 
return  pronounces  acceptance  and  the  algorithm 
procedes  to  the  next  contour  pair. 

Computer  responds 
MAX.  WARP  ANGLE: 

Here  the  user  sets  the  Maximum  Warp  Angle 
described  in  chapter  5.  Default  is  A5  degrees. 

Computer  responds 
START  WITH  WHICH  LEVEL? 

After  which  the  user  responds  with  an  integer, 
n.  The  tektronix  scope  then  displays  the 
contour  line(s)  of  level  n,  labelled  with 
global  loop  numbers.  A carriage  return  will 
cause  level  n+]  to  be  displayed;  a "control-2" 
(end  of  data)  returns  control  to  BRANCH>. 

Again,  the  computer  asks, 

START  WITH  WHICH  LEVEL? 

The  user  responds  with  an  integer, n.  This 
initiates  branching  and  triangulation.  Each 
pair  of  neighboring  contour  levels  is  consi- 
dered respectively  ranging  from  n to  the  last 
level  read  in.  All  contours  from  the  two 
relevant  levels  are  displayed,  and  identified 
by  a number  followed  by  T(top)  or  B (bottom). 

The  user  selects  from  these  loops  the  ones 
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SINGLE  - 


CAP  - 


which  belong  together  for  triangulation. 

(This  need  only  be  used  for  complex  cases  of 
branching.  Normally,  the  AUTOMATIC  command 
handles  branching  adequately.)  The  computer 
asks : 

BOTTOM  LOOP (S)  (LOCAL): 

Here,  the  user  enters  the  loop  numbers  of  all 
the  loops  that  are  to  participate  in  triangu- 
lation, The  loop  numbers  are  delimited  by 
carriage  returns.  The  series  is  terminated  by 
by  a double  carriage  return,  after  which  the 
computer  asks: 

TOP  LOOP (S)  (LOCAL): 

The  same  procedure  is  again  followed  in  input- 
ting the  top  loops.  This  sequence  of  loop 
numbers  is  loaded  into  an  array  called 
Loopstack,  which  is  passed  to  SUBROUTINE  QUAD 
for  triangulation.  If  there  are  more  than  one 
loop  on  either  top  or  bottom  (i.e.  branching) 
QUAD  displays  all  the  loops  and  inquires  of  the 
user  how  to  assemble  them  into  one  loop  by 
asking 

ENTER  LOOP  SEQUENCE: 

The  user  now  enters  the  sequence  in  which  the 
loops  are  to  be  re-constructed  (as  explained 
in  chapter  3).  If  not  all  loops  are  intended 
to  be  included,  the  sequence  may  be  terminated 
with  a control -Z.  Next,  the  computer  asks  how 
to  interconnect  the  loops  by  typing:  (where 
the  loop  sequence  is  i,j,k) 

LiP? 

Enter  the  node  on  loop  i to  be  joined  to  loop  j 
LjP? 

Enter  the  node  on  loop  j to  be  joined  to  loop  i 
etc . 

When  this  series  of  questions  is  answered,  the 
routine  re-arranges  the  branching  loops  into 
one  loop  (as  explained  in  chapter  3),  and 
proceeds  to  the  final  command  prompt  of 
TR IANGULATE>,  Commands  for  that  prompt  are 
discussed  later. 

This  command  performs  a function  similar  to 
that  of  MANUAL.  The  difference  is  that 
global  loop  numbers  are  used  instead  of  local. 
This  enables  any  loop(s)  to  be  defined  as  top 
or  bottom.  The  SINGLE  command  is  used  only 
in  irregular  situations  involving  concentric 
loops.  Execution  procedes  as  in  MANUAL. 

CAP  invokes  an  algorithm  to  form  a pyramidal 
cap  on  a specified  loop.  The  computer  asks: 
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GLOBAL  LOOP  NUMBER: 

To  find  out  what  the  global  loop  number  is 
(in  case  the  user  doesn't  know)  use  the  INSPECT 
command.  This  should  be  done  prior  to 
invoking  the  CAP  command. 

WARP,  CAP,  INSPECT,  and  TOTALS  return  to  the  BRANCH*  prompt. 
AUTOMATIC  returns  to  the  READ>  prompt  after  triangulating  all  avai- 
lable contours. 

The  final  prompt,  TR IANGULATE> , is  only  encountered  during  the 
BRANCH  command,  or  possibly  during  the  AUTOMATI C/POST-ED  IT  command 
(if  triangulation  is  faulty). 

TR I ANGULATE>  responds  to  these  commands: 

AUTOMATIC  - Same  as  the  AUTOMATIC  command  in  the  BRANCH> 
menu . 

DENSITY  - As  default,  the  program  labels  every  other  node 

for  identification.  If  this  labelling  is  too 
dense  to  decipher  (or  too  sparce) , the 
DENSITY  command  provides  help. 

The  computer  asks, 

I ND3=  , where  IND3  is  the  3rd  DO  LOOP 
parameter  in  the  labelling  loop. 

TOTALS  - Same  as  before. 

PART  - The  MOVIE.BYU  graphics  package  enables  panels 

to  be  grouped  together  as  parts.  Each  part 
can  then  be  addressed:  seperately  with  its 
own  set  of  display  parameters.  Here,  the 
PART  command  funnels  all  succeeding  panels  into 
a specified  part  number.  The  computer  asks: 
CHANGE  PART  FROM  n TO  : 

where  n is  the  currently  assigned  part  number. 
Initial  part  use  must  proceed  in  numeric  order. 
That  is,  part  2 must  be  used  before  part  3, 
etc.  However,  once  used,  there  is  no  restric- 
tion on  future  use,  i.e.  the  part  number  could 
be  changed  from  5 to  2,  for  example.  The 
algorithm  safegaurds  these  rules,  and  notifies 
the  user  of  violations, 

UNMAP  - Sets  a flag  that  causes  the  mapping  algorithm 

to  be  bypassed,  so  contours  are  considered  in 
their  initial  orientation  for  triangulation. 
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MAP  - 
ONE  - 

INTERACTIVE  - 


ERASE  - 

RENUMBER  - 


Opposite  of  UNMAP,  Default  is  MAP. 

Permits  the  user  to  manually  define  a single 
panel  from  any  3 (or  k)  nodes  from  the  contour 
pair. 

Permits  interactive  guidance  of  the  triangula- 
tion by  requesting  limits  for  both  loops  between 
which  to  triangulate.  The  computer  types 
TOP:  n~m 

where  n is  the  last  node  to  be  triangulated 
and  m is  the  number  of  nodes  on  the  top  contour, 
(n  and  m are  merely  stated  for  clarity).  The 
user  responds  with  two  integers  delimited  by  a 
comma.  The  first  delimiter  is  generally  n, 
and  the  second  is  always  <_  m.  Delimiters  are 
also  requested  for  the  bottom  loop.  Triangula- 
tion occurs  within  the  prescribed  limits  and 
is  displayed  for  approval.  If  approval  is 
denied,  only  the  newest  panels  are  erased,  and 
the  user  is  invited  to  try  again.  A control-2 
returns  control  to  TRIANGULATE>. 

Erases  al 1 panels  from  the  present  pair  of 
contours. 

Tr  iangul at  ion  is  facilitated  if  node  1 on  the 
top  loop  neighbors  node  1 on  the  bottom  loop. 
RENUMBER  enables  this  by  asking 
NEW  11  NODE 
(TOP  CONTOUR) 

The  renumbered  contours  are  then  displayed. 

This  renumbering  must  occur  before  any  triangu- 
lation. AUTOMATIC  automatically  renumbers 
before  triangulating. 


CONTINUE  - 


Returns  control  to  the  main  program. 
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Dennis  L.  Fung 
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Box  X-  Oak  Ridge,  TN  37830 


Chuck  Weger 
CUNY/UCC 

555  West  57th  Street 
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Department  of  Chemistry 
Simmons  College 
300  Fenway 
Boston,  MA  02115 

fr  McKay  Anderson 

Hercules  Incorporated 
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Albuquerque,  NM  87115 
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Clark  Equipment 
Box  547 

Benton  Harbor,  MI  49022 

Bruce  E.  Koopman 
Sandia  Laboratories 
Division  8323 
Livermore,  CA  94550 
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TRW  DSSG,  Design  Analysis  Center 
One  Space  Park,  M3/2042 
Redondo  Beach,  CA  90278 

^ H W Vannier 

Wenner  Gren  Laboratory 
University  of  Kentucky 
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Ralph  A Kinney 

Electrical  Engineering  Dept. 
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